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ABSTRACT 



Solving an integer least squares (ILS) problem usually consists of two 
stages: reduction and search. This thesis is concerned with the reduction 
process for the ordinary ILS problem and the ellipsoid-constrained ILS prob- 
lem. For the ordinary ILS problem, we dispel common misconceptions on the 
reduction stage in the literature and show what is crucial to the efficiency 
of the search process. The new understanding allows us to design a new re- 
duction algorithm which is more efficient than the well-known LLL reduction 
algorithm. Numerical stability is taken into account in designing the new re- 
duction algorithm. For the ellipsoid-constrained ILS problem, we propose a 
new reduction algorithm which, unlike existing algorithms, uses all the avail- 
able information. Simulation results indicate that new algorithm can greatly 
reduce the computational cost of the search process when the measurement 
noise is large. 
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ABREGE 



La resolution de problemes de moindres carres en nombres entiers (ILS) 
comprend habituellement deux stages: la reduction et la recherche. Cette these 
s'interesse a la reduction pour le probleme ILS ordinaire et le probleme ILS sous 
contrainte d'ellipse. Pour le probleme ILS ordinaire, nous dissipons des erreurs 
communes de comprehension a propos de la reduction dans la litterature et 
nous montrons ce qui est reellement crucial pour l'efficacite de la recherche. Ce 
resultat nous permet de developper un nouvel algorithme de reduction plus effi- 
cace que le celebre algorithme LLL. La stabilite numerique est prise en compte 
dans le developpement du nouvel algorithme. Pour le probleme ILS sous con- 
trainte d'ellipse, nous proposons un nouvel algorithme de reduction qui, con- 
trairement aux algorithmes existants, utilise toute l'information disponible. 
Les resultats de simulations indiquent que le nouvel algorithme reduit con- 
siderablement les couts de calcul de la recherche lorsque la variance du bruit 
est large dans le mo dele lineaire. 
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CHAPTER 1 
Introduction 



Suppose we have the following linear model 

y = Ax + v, (1.1) 

where y £ M m is a measurement vector, A £ M mxri is a design matrix with 
full column rank, x £ IR n is an unknown parameter vector, and v £ M n is a 
noise vector that follows a normal distribution with mean and covariance 
a 2 I. We need to find an estimate x of the unknown vector x. One approach 
is to solve the following real least squares (LS) problem 

min \\y - Ax\\\. (1.2) 



We refer to ( jl.2jl as the standard form of the LS problem. The real least 



squares solution is given by (see, e.g. 



20 



Section 5.3]) 



x = (A T A)- 1 A T y. (1.3) 

Substituting ( 11. ip into ( 11. 30 . we get 

x = x + {A T A)- 1 A T v. (1.4) 

Let W x £ M. nxn be the covariance matrix of x. Using the law of covariance 
propagation on ( II .40 (see, e.g., [33|, p. 329]), it follows that W x = a 2 (A T A)^ 1 . 

In many applications, x is constrained to some discrete integer set T>, 
e.g., a box constraint. Then, one wants to solve the integer least squares (ILS) 
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problem 



mm\\y — Ax\\ 2 . (1-5) 



If T> is the whole integral space Z n , then we refer to 

min \y — Aa;^ (1-6) 



as the ordinary integer least squares ( OILS) problem. In lattice theory, the set 
C(A) = {Ax : x G Z™} is referred to as the lattice generated by A. The OILS 
problem is to find the point in C(A) which is closest to y. For this reason, 
the OILS problem is also called the closest point problem. Since the residual 
y — Ax is orthogonal to the range of A, we have 



\\y - Ax\\ 2 2 = \\y -Ax- A(x - x)\\l 

= \\y- Ax\\l + \\A(x-x)\\ 2 2 

= \\y - Ax\\l + (x — xf A T A(x - x) 

= \\y - Ax\\l + a 2 (x - x) T W^ l (x - x). 

Thus, the OILS problem (11.61) is equivalent to 



;i-7) 



m\n(x — x) T W ^}{x — x). (1. 



We refer to (11.81) as the quadratic form of the OILS problem. 

ILS problems arise from many applications, such as global navigation 
satellite systems, communications, bioinformatics, radar imaging, cryptogra- 
phy, Monte Carlo second- moment estimation, lattice design, etc., see, e.g., [l| 



and 



22) . Unlike the real least squares problem (II. 2p . the ILS problem (II. 5p is 



NP-hard (see [4] and [29(). Since all known algorithms to solve the ILS prob- 
lem have exponential complexity, designing efficient algorithms is crucial for 
real-time applications. A typical approach for solving an ILS problem consists 
of two stages: reduction and search. The reduction process transforms (II. 5p 
to a new ILS problem, where A is reduced to an upper triangular matrix. The 
search process searches for the solution of the new ILS problem in a geometric 



2 



region. The main goal of the reduction process is to make the search process 
more efficient. 

For solving (11.61) . two typical reduction strategies are employed in practice. 
One is the Korkine-Zolotareff (KZ) reduction (see [24]), which transforms the 
original OILS problem to a new one which is optimal for the search process. 
The other is the Lenstra-Lenstra-Lovasz (LLL) reduction (see U), which 
transforms the original OILS problem to a new one which is approximately 
optimal for the search process. Unlike the KZ reduction, it is known how to 
compute the LLL reduction efficiently in polynomial time. For this reason, the 
LLL reduction is more widely used in practice. Simulations in [l| suggest to 
use the KZ reduction only if we have to solve many OILS problems with the 
same generator matrix A. Otherwise, the LLL reduction should be used. For 
the search process, two common search strategies are the Pohst enumeration 



strategy (see 



18) ) and the Schnorr-Euchner enumeration strategy (see 



Both examine the lattice points lying inside a hyper-sphere, but in a different 
order. Simulations in [lj] indicate that the Schnorr-Euchner strategy is more 
efficient than the Pohst strategy. 

In high precision relative global navigation satellite systems (GNSS) posi- 
tioning, a key component is to resolve the unknown so-called double differenced 
cycle ambiguities of the carrier phase data as integers. The most success- 
ful method of ambiguity resolution in the GNSS literature is the well-known 
LAMBDA (Least-squares AMBiguity Decorrelation Adjustment) method pre- 



sented by Teunissen (see, e.g., 



34 



35 



36 



38 



39) ). This method solves the 



OILS problem (11.81) . A detailed description of the LAMBDA algorithm and 
implementation is given by [16) . Its software (Fortran version and MATLAB 
version) is available from Delft University of Technology. Frequently asked 
questions and misunderstanding about the LAMBDA method arc addressed 
by 171 ] . Recently, a modified method called MLAMBDA was proposed by 
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Chang etc in [h]], which was then further modified and extended to han- 
dle mixed ILS problems by using orthogonal transformations, resulting in the 



MATLAB package MILS (see [13]) 



In some applications, the point Ax in (I1.5P is constrained to be inside a 
given hyper-sphere. The ILS problem becomes 



min||y- Ax\\%, £ = {x e Z n : \\Ax\\ 2 2 < a 2 }, (1.9) 



see |9j and 15]. We refer to ( II. 9p as the ellipsoid-constrained integer least 
squares (EILS) problem. To solve the EILS problem, [9J proposed the LLL 
reduction in the reduction stage and modified the Schnorr-Euchner search 
strategy to handle the ellipsoidal constraint in the search stage. 

The contribution of this thesis is two-fold. The first is to dispel com- 
mon misconceptions on the reduction process appearing in the ILS literature. 
The second is to present more efficient algorithms to solve the ILS and EILS 
problem. The thesis is organized as follows. 

In Chapter [2J we review the typical methods to solve an OILS problem. 
Specifically, we introduce the LLL reduction method and the Schnorr-Euchner 
search strategy. 

In Chapter [3j we discuss the typical methods to solve the quadratic form 
of the OILS problem. Our focus will be on LAMBDA reduction and the 
modified reduction (MREDUCTION) given in [lo|. 

According to the literature, there are two goals the reduction process 
should achieve to make the search process efficient. One of them is to trans- 
form the generator matrix A in (11. 6p or the covariance matrix W x in ( II. 8p to 
a matrix which is as close to diagonal as possible. A covariance matrix which 
is close to diagonal means that there is little correlation between its random 
variables. In other words, the goal of the reduction in the GNSS context is 
to decorrelate the ambiguities as much as possible. To achieve this goal, the 
reduction process uses so-called integer Gauss transformations. In Chapter HJ 
we show that, contrary to common belief, this goal will not make the search 
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more efficient. We provide a new explanation on the role of integer Gauss 
transformations in the reduction process. This new understanding results in 
modifications to the existing reduction methods. Numerical simulations indi- 
cate that our new algorithms are more efficient than the existing algorithms. 
Finally, we discuss another misconception in some GNSS literature where it 
is believed that the reduction process should reduce the condition number of 
the covariance matrix. We provide examples that show that this goal should 
be discarded. 

In Chapter 0, we discuss the EILS problem. One drawback with the 
existing reduction algorithms is that the search time becomes more and more 
prohibitive as the noise v in (11. ip gets larger. We present a new reduction 
algorithm which, unlike the LLL reduction, uses the information of the input 
vector y and the ellipsoidal constraint. Then, we provide simulations that 
show that our proposed approach greatly improves the search process for large 
noise. 

Finally, we summarize our results and mention our future work in Chapter 

El 

We now describe the notation used in this thesis. The sets of all real 
and integer m x n matrices are denoted by R mxn and Z mxn , respectively, and 
the set of real and integer n- vectors are denoted by W 1 and Z n , respectively. 
Bold upper case letters and bold lower case letters denote matrices and vectors, 
respectively. The identity matrix is denoted by J and its ith column is denoted 
by ej. Superscript T denotes the transpose. The 2-norm of a vector or a 
matrix is denoted by || ■ 1 1 2 . MATLAB notation is used to denote a sub matrix. 
Specifically, if A = (a^) G M mxn , then A(i, :) denotes the zth row, A(:, j) the 
jth column, and A(ii : 2*2, jfi : J2) the submatrix formed by rows i\ to 12 and 
columns ji to j%. For the element of A, we denote it by or A(i,j). 
For the ith entry of a vector a, we denote it by a«. For a scalar z G R, we 
use \_z] to denote its nearest integer. If there is a tie, \_z] denotes the one 
with smaller magnitude. The operation sign(z) returns — 1 if z < and 1 if 
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z > 0. For a random vector x G W l , x ~ A/"(0, <r 2 /) means that as follows a 
normal distribution with mean and covariance matrix a 2 I. We use i.i.d. to 
abbreviate "independently and identically distributed". 
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CHAPTER 2 
The OILS problem in the standard form 



The OILS problem is defined as 

mjn||y-Ax||*, (2.1) 

where y G M n and A G M mxn has full column rank. In this chapter, we 
review the typical methods to solve the OILS problem. A typical approach for 
solving an OILS problem consists of two stages: reduction and search. The 
main goal of the reduction process is to make the search process efficient. In 
order to better understand the aims of the reduction, we first introduce the 
search process and the Schnorr-Euchner search strategy in Section 12.11 Then, 
we present the reduction process and the well-known LLL reduction in Section 

m 

2.1 Search process 

Suppose that after the reduction stage, the OILS problem (12.11) is transformed 
to 

min \\y- Rz\\l, (2.2) 
zez"" y ll2 ' V 1 

where R G ~R nxn is nonsingular upper triangular and y G M. n . Assume that 
the solution of ( 12. 21) satisfies the bound 



\y-Rz\\l<{3 
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or equivalent ly 

n n 

^{yk - ^2 Tk i z i ~ r kkZkf < P 2 - (2.3) 

k=l j=k+l 

Since ( 12. 3 p is a hyper-ellipsoid, we refer to it as a search ellipsoid. The 
goal of the search process is to find an integer point inside the search ellipsoid 
which minimizes the left-hand side of (12. 3p . 

Let 

n 

c n = yn/r n n, c k = (y k - r kjZj) /r kk , k = n - (2.4) 

j=k+l 

Note that c k depends on z k+1 , . . . , z n . Substituting ( 12. 4 p in (12. 3p . we have 

J2r 2 kk (zk - c k f < (3\ (2.5) 

k=l 

If z satisfies the bound, then it must also satisfy inequalities 



level n : r 2 nn (z n - c n f < (3 2 



(2.6) 



level k : r 2 kk (z k - c k f < (3 2 - ^ 

i=k+l 



(2.7) 



levellrr^^i-ci)^^ 2 -^- 



i=2 



(2.8) 



The search process starts at level n and moves down to level 1. At level 
k, z k is determined for k = n : — 1 : 1. From (12.70 . the range of z k is [l k ,u k ], 
where 



and 



u k 



c k -(/3 2 - rate-*, 

i=k+l 



2U/2 



/Vkk\ 



Ck + (P 2 - r*(*-*) 2 ) 1/2 /\ru.\ 

i=k+l 



8 



There are two typical strategies to examine the integers inside [lk,Uk\- In 
the Pohst strategy (see If]]), the integers are chosen in the ascending order 

h, h + 1, h + 2, . . . , Uk- 



32]), the integers are chosen in the zig- 



In the Schnorr-Euchner strategy (see 
zag order 

{[c k ] , [c k ] - 1, [c k ] + 1, [c k ] - 2, ... , if c fc < [c k ] / 
(2.9) 
[c k ~\ , [c k ~\ + 1, [c k ] - 1, [c k ] + 2, . . . , if c fc > [c k ] ■ 

Observe that in Schnorr-Euchner strategy, once an integer Zk does not 
satisfy (12.71) . all the following integers in the sequence will not satisfy it. These 
integers can be pruned from the search process. Such a property does not exist 
in the Pohst strategy. Another benefit with the Schnorr-Euchner enumeration 
order is that the first points examined are more likely to minimize (12.51) than 
the last points examined. As will be seen in the next paragraph, this allows to 
shrink the search ellipsoid faster. Simulations in 1[ confirm that the Schnorr- 
Euchner strategy is more efficient than the Pohst strategy. 

We now describe the search process using the Schnorr-Euchner strategy. 
At level n, we compute c n by (12.4ft and set z n = |_c n ] . If (12. 6ft is not satisfied, 
no integer can satisfy (12. 5p . Otherwise, we go to level n — 1, compute c n _i and 
set 2„_i = L c n-i]- If (12.71) does not hold, we go back to level n and choose 
z n to be the second nearest integer to c n . Otherwise, we move down to level 
n — 2. When we reach level 1, we compute C\ and set z\ = |_Ci] . Then, if (12. 8p 
is satisfied, we set z — [zi, . . . , z n ] T , where z is a full integer point inside the 
search ellipsoid. We update f3 by setting (3 2 = Ylk=i r kk(^k ~ c k) 2 - This step 
allows to eliminate more points by "shrinking" the search ellipsoid. Now we 
search for a better point than z. If one is found, we update z. We move up 
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k = 3 



k = 2 



k = 1 



Figure 2-1: A three-dimensional example of a search tree 

to level 2 and choose ^2 to be the next nearest integer to C2, where "next" is 
relative to z-i- If inequality ( 12.7ft holds at level 2, we move down to level 1 and 
update Z\\ otherwise, we move up to level 3 and update z 3 . The procedure 
continues until we reach level n and ( 12.61) is not satisfied. The last full integer 
point found is the OILS solution. The described search process is a depth-first 
search, see an example of a search tree in Fig. 12-11 The root node does not 
correspond to an element in z, but is used to unite the branches into a tree. 
Note that the integers are enumerated according to order (I2.9p when we move 
up one level in the search. If the initial value of /3 is oo, the first integer 
point found is called the Babai integer point. We describe the process as an 
algorithm, see |7|. 

Algorithm 2.1.1. (SEARCH) Given nonsingular upper triangular matrix 
R g R nxn and y G M n . The Schnorr-Euchner search algorithm finds the 
optimal solution to min z6 zn \\y — 
function: z = SEARCH (iZ, y) 

1. (Initialization) Set k = n, (3 = oo. 

2. Compute Ck from (I2.4H . Set Zk = [ck], = sgn(cfc — Zk) 

3. (Main step) 
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go to Step 4 
else if k > 1 

k — k — 1, go to Step 2 
else / / case k = 1 

go to Step 5 
end 

4. (Invalid point) 
if A; = n 

terminate 
else 

k = k + 1, go to Step 6 
end 

5. (Found valid point) 

Set z = z, P 2 = E n k=1 r 2 kk (z k -c k y 
k = k + 1, go to Step 6 

6. (Enumeration at level k) 

Set z k = z k + A k , A k = -A fc - sgn(A fc ) 
go to Step 3. 

2.2 Reduction process 

In the reduction process, we transform A to an upper triangular matrix R 
which has properties that make the search process more efficient. Since (12. ip 
uses the 2-norm, we can apply an orthogonal transformation Q T to the left 
of A as long as it is also applied to the left of y, i.e., \\y — Ax\\?, = \\Q T y — 
Q T Ax\\2- In order to maintain the integer nature of x, the transformation 
Z applied to the right of A must be an integer matrix, whose inverse is also 
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an integer matrix. It is easy to verify that |det(-Z)| = 1, as both det(Z') and 
det(Z _1 ) are integers, and det(Z)det(Z' _1 ) = 1. Such integer matrices are 
referred to as unimodular matrices. 

The transformations on A can be described as a QRZ factorization of A: 



Q T AZ 



R 




or A = QTRZ \ (2.10) 



where Q = [Q 1 ,Q 2 ] £ M. mxm is orthogonal, R G M. n is nonsingular upper 
triangular and Z G Z nxn is unimodular (see 9|). Using this factorization, 

\\y - Ax\\l = \\Q\y - RZ~ x x\& + WQlvWl (2.11) 

Let 

y = Q T lV} z = Z- l x. (2.12) 
Then, we can rewrite (12. ip as 

min||y-.R*||!. (2.13) 

If z is the solution of the transformed OILS problem ( 12. 13ft . then x = Zz is 
the solution of the original OILS problem ( 12. ip . 

Note that R in (12.101) is not unique. A different Z usually leads to a 
different R. Without loss of generality, we assume that the diagonal entries 
of R are positive in this thesis. Certain properties of R can make the search 
process much more efficient. If R is diagonal, then simply rounding all Ck to 



the nearest integer gives the optimal solution. In 



30], it is claimed that as R 



gets closer to a diagonal matrix, the complexity of the search process decreases. 



In Section 14.11 we show that this claim is not true. The crucial property that 
R must strive for is that should be as large as possible for large k. We 
motivate this property by a two dimensional case (n = 2). Let r 2 2 <C r\x- 
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In the search process, the bound ( 12. 6 j) at level 2 is very loose, implying that 
there are many valid integers z^- However, the bound f |2.8j) at level 1 is very 
tight. Hence, after Z2 is fixed, there is a high probability that no valid integer 



Z\ exists. We must enumerate many integers at level 2 before we can fine 
valid integer at level 1. This is the so-called search halting problem (see 



a 



Note that det(A T A) = det(R T R) = r\ x ...r 2 nn is constant, independent of 
the choice of Z. Making r kk as large as possible for large k implies making 
as small as possible for small k. Hence, we can say that the reduction process 
should strive for 

ru«...«r nn . (2.14) 
The typical reduction method for the OILS problem is the LLL reduction 



see 



251 ] and 22J). The LLL reduction can be written in the form of a QRZ 
factorization, where R satisfies the following criteria 

< 2 r k-i.fc-i» J = k,---,n (2.15) 
r fc -i,*-i < 6y/ri_ ltk + r 2 kk , 1<5<2 and k = 2,...,n. (2.16) 

Notice that taking 5 = 1 is best to strive for ( I2.14p . In this thesis, we always 
take 5 = 1. The LLL reduction cannot guarantee r n < ... < r nn , but 
substituting ( I2.15P into (I2.16p . we see that it can guarantee 

2 

5"k-i,fc-i < — f=rkk, k = 2,...,n. 

In our implementation of the LLL reduction, we use two types of unimod- 
ular transformations to ensure properties (I2.15P and (I2.16p . These are integer 
Gauss transformations and permutation matrices. In the following, we present 
the effect of these transformations on the R factor of the QR factorization of 
A. 
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2.2.1 Integer Gauss transformations 



Suppose we are given an upper triangular matrix R with positive diagonal 
entries. An integer Gauss transformation (IGT) Z^ has the following form 

Zij = I — /ie^ej, jj, is an integer. (2-17) 

Applying (i < j) to R from the right gives 

R = RZ l3 = R - fiRe.ej. 

Thus R is the same as R, except that 

Taking \i = [r^/ru] ensures that |fy| < \ra. Similarly, an IGT [i > j) 
can be applied to a unit lower triangular matrix L from the right to ensure 
that \kj\ < 1/2. 

2.2.2 Permutations 

For R to satisfy (12.161) with 5 = 1, we sometimes need to permute its columns. 



In the reduction process, if r^-i^-i > yr^_ lfe + r| fc , we permute columns k 
and k — 1. 



RPk~i, 



R\i Rn Ris 
R22 R23 
R33 



where 



-Pfc-i, 



- k-2 



■ n—k 






1 






^fc-l,fc-l 






, R22 — 




1 







r kk 
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R 12 = [R(l :k — 2,k — l), R(l :k-2, k)}. 

As a result, R is no longer upper triangular. To make it upper triangular 
again, we apply a Givens rotation G to zero element in R 2 2- 



GR22 — R22, o r 



c s 
-s c 



Tk-l,k fk-l,k-l 

r k k 



Tk-l,k-l r k-l,k 

ffc fc 



where 



f k -l,k-l 



+ 7tfc, c 



k-l,k ~ ' kki 



Tk-l,k 
Tk-hk-1 



r k k 



Tk-1M-1 



r fc _i ifc = crfc_i,fc_i, r fcfc = -srfe_i,fc_i. 



(2.18) 
(2.19) 



Therefore, we have 



Qk-l,kRPk-l,k — R 



Rn R12 R13 
R22 R23 
R33 

R23 — GR 



■k-2 



G 



■ n—k 



23- 



After the permutation, inequality Tk-i,k-\ < y r k-i,k + r W now holds. While 
a permutation does not guarantee fk-\,k-\ < r kk , it does guarantee that 

fk-x,k-\ < r k-i,k-i and f k k > r kk . 

Such a permutation is useful since the diagonal elements of R are now closer 
to (EH}. 

2.2.3 LLL reduction 



The LLL reduction algorithm given in 



9J starts by finding the QR decompo- 



sition of A by Householder transformations, then computes y and works with 
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R from left to right. At the kth column of R, the algorithm applies IGTs to 
ensure that |r ifc | < \ra ior i = k — 1 : — 1:1. Then, if inequality (I2.16P holds, 
it moves to column k + 1; otherwise it permutes columns k and k — 1, applies 
a Givens rotation to R from the left to bring R back to an upper triangular 
form, simultaneously applies the same Givens rotation to y, and moves to col- 
umn k — 1. We describe our implementation of the LLL reduction as follows 
(see jsl). 

Algorithm 2.2.1. (LLL Reduction). Given the generator matrix A G K mxn 
and the input vector y G M m . The algorithm returns the reduced upper 
triangular matrix R G R nxn , the unimodular matrix Z G Z nxn , and the 
vector y G M n . 

function: [J2, Z, y] = LLL(A, y) 

Compute QR factorization of A and set y = Qjy 

Z = I 

k = 2 

while k < n 

for z = — 1 : — 1 : 1 

Apply IGT Z ik to R, .i.e., R = RZ lk 
Update Z, i.e., Z = ZZ^ 

end 



if r fc _i, fe _i > yjri_ 1>k + rjk 

Interchange columns k and — 1 of R and Z 

Transform R to an upper triangular matrix by a Givens rotation 

Apply the same Givens rotation to y 

if k > 2 

k = k-l 
end 
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else 

k = k + 1 

end 

end 
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CHAPTER 3 
The OILS problem in the quadratic form 



A prerequisite for high precision relative GNSS positioning is to resolve the 
unknown double differenced cycle ambiguities of the carrier phase data as 
integers. This turns out to be an OILS problem. Suppose i G 1" is the 
real- valued least squares estimate of the integer parameter vector igZ™ (i.e., 
the double differenced integer ambiguity vector in the GNSS context) and 
W x € R nxn is its covariance matrix, which is symmetric positive definite. 
The OILS estimate x is the solution of the minimization problem: 



mm 



m(x — x) T W \ l {x — x). 



(3.1) 



Although (13.11) is in the form of an integer quadratic optimization prob- 
lem, it is easy to rewrite it in the standard OILS form ( 12. II) . We refer to 
(13. ip as the quadratic form of the OILS problem. In Section I3.1[ we discuss 



the reduction process used in the GNSS literature. In Section 13.21 we review 



the reduction stage in the LAMBDA 



Adjustment) method (e.g., 



34 



35 



36 



east-squares AMBiguity Decorrelation 



38 



39[). In Section [3731 we introduce 



the improvements to the reduction provided by the MLAMBDA (Modified 



LAMBDA) method (see 



10] ). Finally in Section 13. 4[ we briefly review the 



search process in the quadratic form of the OILS problem. 
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3.1 Reduction process 

The reduction step uses a unimodular matrix Z to transform ( 13. II) into 

min(2 - z) T W^(z - z), (3.2) 

where z = Z T x, z = Z T x and W z = Z T WxZ. If z is the integer minimizer 
of (13.21) . then x = Z z is the integer minimizer of ( 13. ip . The benefit of the 
reduction step is that the search in the new optimization problem ( 13 .2[) can 
be much more efficient. If W z is a diagonal matrix, then the transformed 
ambiguities Zi,...,z n are uncorrelated to each other. In this case, simply 
setting Zi = [zfl , for i = 1 : n, would minimize the objective function. 
Let the L T DL factorization of W z be 

W z = L T DL, (3.3) 

where L is unit lower triangular and D = diag(o?i, . . . , d n ) with di > 0. These 
factors have a statistical interpretation. Let Z{ denote the least-squares es- 
timate of Zi when z i+ i, . . . , z n are fixed. As shown in 38|, p. 337], di is the 
variance of Zi, which is denoted by cr|.. Furthermore, Uj = (JziZjOHf for i > j, 
where a ZiZi denotes the covariance between zt and Zj. 



In the literature (see, e.g., 



16|] 



33j, p. 498] and [3J, p. 369]), it is often 



mentionned that the following two goals should be pursued in the reduction 
process because they are crucial for the efficiency of the search process: 
(i) W z is as diagonal as possible. From ( 13.31) . for i 7^ j, making L(i + 1 : n, 1) 
and L(j + 1 :n, j) closer to makes W z (i, j) closer to 0. Hence, making 
the absolute values of the off-diagonal entries of L as small as possible 
makes W z as diagonal as possible. A covariance matrix which is close 
to diagonal means that there is little correlation between its random 
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variables. In other words, the goal of the reduction is to decorrelate the 
ambiguities as much as possible, 
(ii) The diagonal entries of D are distributed in decreasing order if possible, 
i.e., one strives for 

di > d 2 > ■ • • > d n . (3.4) 

Note that we want d\ to be as large as possible and d n to be as small 
as possible. We can show that striving for (13. 4p is equivalent to striving for 
(12.141) . where di corresponds to r^ 2 for i = 1 :n. In the reduction process of the 
LAMBDA method, the unimodular matrix Z is constructed by a sequence of 
integer Gauss transformations and permutations. The reduction process starts 
with the L T DL factorization of W& and updates the factors to give the L T DL 
factorization of W ~ z . The main contribution of this thesis is to show that, 
contrary to common belief, the first goal will not make the search process 
more efficient. While lower triangular integer Gauss transformations are used 
to make the absolute values of the off-diagonal entries of L as small as possible, 
we argue that they are useful because they help achieve the second goal. 



In 



261 ]. [271 ] and [431 ]. instead of (i) and (ii), the condition number of W 



is used to evaluate the reduction process. In Section I4.6[ we show that this 
criterion can be misleading and that it is not as effective as (ii). 

3.1.1 Integer Gauss transformations 

Integer Gauss transformations (IGTs) were first introduced in Section 12.2. 1[ 
We apply Z^ with /i = [1^1 (see (I2.17P ) to L from the right, i.e., L = LZij, 
to make as small as possible. This ensures that 

|y<l/2, i>j. (3.5) 
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We use the fo 



lowing algorithm to apply the IGT Zj,- to transform the OILS 



problem (see [lOl]). 

Algorithm 3.1.1. (Integer Gauss Transformations). Given a unit lower tri- 
angular L G IR nxn , index pair x E M. n and Z E Z nxn . This algorithm 
applies the integer Gauss transformation to L such that | (LZ) (i, j)\ < 1/2, 
then computes Z^x and ZZy, which overwrite x and Z, respectively, 
function: [L,x,Z] = GAUSS (L, i, j, x, Z) 

if p. ^ 

L(i : n, j) = L(i : n, j) - fiL(i : n, i) 
Z(l:n,j) = Z(l:n,i)-^(l:n,i) 

end 



3.1.2 Permutations 

In order to strive for order (13.41) . symmetric permutations of the covariance 
matrix are needed in the reduction. After a permutation, the factors L 
and D of the L T DL factorization have to be updated. 

If we partition the L T DL factorization of as follows 



L 1 DL 



-^11 ^21 lj 



-^22 "^32 



'33 



D 1 



D ? 









k- 


1 


L21 


L22 




2 




L31 


L32 


£33 


n — 


fc-1 




2 


n—k— 


1 





Let 



1 

1 



• fc-i 



fc.fc+i 



■n-fc-l 
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It can be shown that P\ k+1 W k,k+i has the L T DL factorization 



Pk,k+l W xPk,k+l 



-^11 



^21 
^22 



-^31 
-^32 



J 33 



D 1 



L31 L32 £33 



(3.6) 



where 



di 



d 



k+l 



d 



k+l 



d k + I 



k+l,kdk+l, 



d. 



fe+i 



'22 



k+l,k 



I 



L2I 
£32 — L32P 



-I 



k+l,k 



k+l,k 



'21 



dk+il 



k+i,k 



d 



k+l 



h+i, 



I 



j t-k+l.k 

d k+l 



(3.7) 



(3.8) 



L(k:k + l,l:k-l), (3.9) 



L(k + 2:n,k + 1) L(k + 2:n,l:k) 



(3.10) 



We refer to such an operation as a permutation between pair (k, k + l). We 



101). 



describe the process as an algorithm (see 
Algorithm 3.1.2. (Permutations). Given the L and D factors of the L T DL 
factorization of W x G IR nxn , index k, scalar 5 which is equal to c4 +1 in (13.71) . 
x G M™, and Z G Z nxn . This algorithm computes the updated L and D 
factors in (13.61) after W&s kth row and (k + l)th row, and kth column and 
(k + l)th column are interchanged, respectively. It also interchanges je's kth 
entry and (k + l)th entry and Z's kth column and (k + l)th column. 

function: [L, D, x, Z] = PERMUTE(L, D, k, 6, x, Z) 

r] = D(k,k)/5 II see fl3~7D 

A = D(k + l,k + l)L(k + 1, k)/5 II see 

D{k,k) = r]D(k + l,k + l) II see fl3T7D 
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L(k:k + l,l:k-l) // see § 



D(k + l,k + l) = 5 

-L(k+l,k) 1 

L(Jfc:Jfe + l,l:Jfe-l) = 

i] X 
L(k + l,k) = X 

swap columns L(k + 2:n,k) and L(k + 2 : n, k + 1) // see (13.1 Op 
swap columns Z(l :n, k) and Z(l :n, k + 1) 
swap entries £e(/c) and x(k + 1) 

3.2 LAMBDA reduction 
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We now describe the reduction process of the LAMBDA method (see 
Sect. 3]). First, it computes the L T DL factorization of W x . The algorithm 
starts with column n — 1 of L. At column k, IGTs are applied to ensure 
that the absolute values of the entries below the (k, k)th entry are as small 
as possible. Then, if > c4+i holds (see (13. 7p ). it moves to column k — 1; 
otherwise it permutes pair (k, k + 1) and moves back to the initial position 
k = n — 1. The algorithm uses a variable (kl in Algorithm 13.2.11 below) to 
track down the columns whose off-diagonal entries in magnitude are already 
bounded above by 1/2 due to previous integer Gauss transformations. We 



present the implementation given in 



io|. 



Algorithm 3.2.1. (LAMBDA REDUCTION). Given the covariance matrix 
W x and real- valued LS estimate x of x. This algorithm computes an integer 
unimodular matrix Z and the L T DL factorization W * z = Z T WxZ = L T DL, 
where L and D are updated from the factors of the L T DL factorization of 
W x . This algorithm also computes z = Z T x, which overwrites x. 
function: [Z, L, D, x] = REDUCTION (W x , x) 
Compute the L T DL factorization of W x : W x = L T DL 
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Z = I 

k = n — 1 
kl = k 
while k > 

if k < kl 

for i = k + 1 : n 

[L, x, Z] = GAUSS(X, i, k, x, Z) 

end 
end 

D(k + 1, A; + 1) = D(k, k) + L(k + 1, k) 2 D(k + 1, k + 1) 
if D(jfe + 1, jfe + 1) < D(k + l,k + l) 

[L, D, x, Z] = PERMUTE(£, D, k, D(k + l,k + 1), x, Z) 

kl = k 

k = n — 1 

else 

k = k-l 

end 

end 

When the reduction process is finished, we have 

iy<l/2, J = 1 A; - 1, (3.11) 

4+i < d k + l 2 k+1)k d k+1 , k = 1,2, . . . ,n - 1. (3.12) 

Note that (13. lip and (13. 12)) are the lower-triangular equivalent of the LLL 
reduction properties (12. 1 5)) and (12.161) with 5 = 1, respectively. 
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3.3 MLAMBDA reduction 



In the MLAMBDA method, several strategies are proposed to reduce the com- 
putational complexity of the reduction process. 



3.3.1 Symmetric pivoting strategy 



In striving for fj3.4j) . LAMBDA reduction performs symmetric permutations 
of the covariance matrix W x . After each permutation, an IGT is needed to 
update the L and D factors of the L T DL factorization of W x . One idea 
is to apply permutations to W x before computing its L T DL factorization. 
This way, the cost of the IGT associated with a permutation is saved. Once 
the L T DL factorization is computed, new permutations are usually needed 
to strive for (13. 4p . Nevertheless, this strategy usually reduces the number of 
permutations done after we have the L and D factors. First, we show how to 
compute the L T DL of W x without pivoting. We partition the W x = L T DL 
as follows 







~ T 










W x q 




L I 




D 




L 


Q T Qnn 




1 








f 1 



We can see that 



dn QVm; ^ Q/ ^"ni a 



ldJ T = L T DL. 



(3.13) 



We recurse on W x — ld n l T to find the complete factorization. Now we intro- 
duce the symmetric pivoting strategy. Since we strive for (13 .4p . we first sym- 
metrically permute the smallest diagonal entry of W x to position (n, n). With 
( 13.131) . we compute d n and I. We continue this procedure with W x — ld n l T . Fi- 
nally, we get the L T DL factorization of a permuted W x . The implementation 



of this strategy is described in 10] 
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Algorithm 3.3.1. (L T DL factorization with symmetric pivoting). Suppose 
W x £ M nxn is symmetric positive definite. This algorithm computes a per- 
mutation P, a unit lower triangular matrix L and a diagonal D such that 
P T W X P = L T DL. The strict lower triangular part of W x is overwritten by 
that of L and the diagonal part of W x is overwritten by that of D. 

P = In 

for k — n : — 1 : 1 

q = argmini<j-<fc W x (j,j) 
swap P(:,k) and P(-,q) 
swap W x (k, :) and Wj(g, :) 
swap W x (:,k) and Wa(:,g) 
W x (k, 1 : fc - 1) = W x (k, 1 : k - 1)/W x (k, k) 
W 4 (l:jfe-l,l:ife-l) = Wa(l:Jfe-l,l:Jfe-l) - W^(A;, 1 :k — l) T * 

W x (k,k)*W x (k,l:k-l) 

end 

3.3.2 Greedy selection strategy 

The reduction process starts with the L T DL factorization with pivoting. In 
order to further reduce the number of permutations, a greedy selection strategy 
is proposed. As shown in Section 13.21 the reduction process of the LAMBDA 
method permute pairs (k, k + 1) from right to left. If for some index k, we have 
dk+i d k and c^+i < d k+ i, then we permute pair (k, k+1) and we move to 
column k+1 of L. Now, it is likely that we also have to permute pair (k + 1, k + 
2), and so on. As a result, it is possible that some of the permutations done 
before reaching index k are wasted. To avoid these unnecessary permutations, 
instead of looping k from n — 1 to 1 as in Algorithm 13.2.11 we choose the 
index k such that (4+i decreases most after a permutation for pair (k,k + 1) 
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is performed. In other words, we first permute the pairs (k, k + 1) for which 
we are most confident of the order. We define k by 



k = arg min {dj + i/dj + i : oL+i < cL+i}. (3-14) 

l<j<n— 1 

If no k can be found, no more permutations are applied. 



3.3.3 Lazy transformation strategy 

In LAMBDA reduction, IGTs can be applied to the same entries in L numerous 
times. We now explain how this can occur. When we permute pair (k, k + 1), 
the entries of L{k : k + 1, 1 : k — 1) are modified (see ( 13 .8 j) ) . If the absolute 
values of the entries of L(k : k + 1, 1 : k — 1) are bounded above by 1/2 before 
the permutation, then these bounds may no longer hold after the permutation. 
Hence, new IGTs have to be applied. To avoid this extra work, we want to 
defer as much IGTs as possible to the end of the reduction process. From 
(13. 7p . if d k < d k+ i, we need |/fc+i,*l to be as small as possible to determine the 



order of the ambiguities. Therefore, at first, we apply IGTs only on some of 
the subdiagonal entries of L. Then, when no more permutations occur, IGTs 
are applied to all the entries in the strictly lower triangular part of L. This 



strategy is called a "lazy" transformation strategy in 
3.3.4 The reduction algorithm 



io| 



We present the modified reduction algorithm (MREDUCTION) given in 
It uses an {n + l)-dimensional vector ChangeFlag to track if h+i,k is modified 
by the last permutation. 

Algorithm 3.3.2. (MREDUCTION) Given the covariance matrix W & and 
real-valued LS estimate x of x. This algorithm computes an integer unimod- 
ular matrix Z and the L T DL factorization Wt- = Z T W^Z = L T DL, where 
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L and D are updated from the factors of the L T DL factorization of Wj. This 
algorithm also computes z = Z T x, which overwrites x. 

function: [Z, L, D, x] = MREDUCTION(V^, x) 

Compute the L T DL factorization of W& 

with symmetric pivoting P T WxP = L T DL 

x = P T x 

Z = P 

Set all elements of ChangeFlag(l:n+l) to ones 
while true 

minratio = 1 

for k — 1 : n — 1 
if D ^ k ) < i 

£>(fc+l,fc+l) 

if ChangeFlag(A; + 1) = 1 

[L, x, Z] = GAUSS(L, k + l,k, x, Z) 

D(k + 1, k + 1) = D(k, k) + L(k + 1, k) 2 D(k + 1, k + 1) 

ChangeFlag(A; + 1) = 1 
end 

, _ r>(fc+i,fc+i) 
tmp — £,( fe+1)fc+1 ) 

if tmp < minratio 

i — k II see fl3T4|) 

minratio = tmp 

d = D(k + 1, k + 1) 
end 
end 

end 

if minratio = 1 

break while loop 
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end 

[L, D, x, Z] = PERMUTE(L, D, i, d, x, Z) 
Set ChangeFlag(i : i + 2) to ones 

end 

/ / Apply IGTs to L's strictly lower triangular part 
for k = 1 : n — 1 

for % = k + 1 : n 

[L, x, Z] = GAUSS(L, i, k, x, Z) 

end 

end 

Numerical simulations in 



10j show that MLAMBDA can be much faster 
than LAMBDA implemented in Delft's LAMBDA package (MATLAB, version 
2.0) for high dimensional problems. In Section fl~5l our numerical simulations 
indicate that MREDUCTION can be numerically unstable on some problems. 
In these problems, MLAMBDA finds a worse solution to the OILS problem 
than LAMBDA. 

3.4 Search process 

After the reduction process, the search process starts. The critical step in 
understanding the search process is to rewrite the objective function (13. 2 p in 
terms of a sum-of-squares, similar to ( 12. 51) . Substituting the L T DL factoriza- 
tion in (I3.2p . we get 

mm(z-z) T L~ l D- l L- T (z-z). (3.15) 

Define z as 

z = z- L- T (z- z), (3.16) 
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or equivalent ly 

L T (z — z) = z — z, 

which can be expanded to 

n 

*i = % + li i( Zi ~ j = n : -I : 1. (3-17) 

1=3+1 

Observe that Zj depends on zj + i,...,z n . With ( 13.1 6p . we can rewrite the 
optimization problem (I3.15P as follows 



min(z - zfD'Hz - z), (3.18) 



or equivalently 



minV 1 ^ Zj) \ (3.19) 



j=l 3 



Assume that the solution of (13.191) satisfies the bound 



3=1 3 

Note that (13.201) is a hyper-ellipsoid, which we refer to as an ambiguity search 
space. If z satisfies (13.201) . then it must also satisfy inequalities 



level n : < p , 



level k . ^ ~ *» <f-J2 {Zl - ^ (3.21) 

i=k+l 1 





-z n f 




d n 






{z k 


-z k f 




d k 






(zi 


~ zi? 



level 1: < (3 2 - ^ 



n i - \2 

Zi Zj) 



i=2 % 
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The search process starts at level n and moves down to level 1. From f )3.2ip . 
the range of is where 



Zk 



=fc+i 



and 



Uk 



1/2 



i=k+l 



(3.22) 



(3.23) 



With the inequalities at each level, the search for the OILS solution can be 
done with the same procedure shown in Section 12.11 
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CHAPTER 4 

Reduction Misconceptions and New Reduction Algorithms 



Basically there are two communities studying ILS problems: the information 
theory and communications community and the GNSS community. Typically, 
the former uses the OILS problem in the standard form, while the later uses 
the quadratic form. In Section [2j we presented the OILS problem in the 
standard form and the LLL reduction method. In Section [31 we presented 
the OILS problem in the quadratic form and the LAMBDA reduction and 
the MREDUCTION methods. It appears that there are two misconceptions 
about the reduction process in the literature. The first is that the reduc- 
tion process should decorrelate the covariance matrix of the real least squares 



estimate as far as possible, i.e., make the off-diagonal entries o 



p. 498] and 



t 
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he covari- 
p. 369]). 



ance matrix as small as possible (see, e.g., [la], [33 
This misconception also appears in the communications literature, where it is 
claimed that the search process will be faster if the reduction process makes 
the off-diagonal entries of the triangular matrix R as small as possible (see, 



e.g. 



9] and [30[). The second is that the reduction process should reduce 



he condition number of the the covariance matrix (see, e.g., 



26] 



27j and 



43|). In this Chapter, we show that both are incorrect in Sections 14. II and | 4.6[ 



respectively. Our results will provide insight on the role of lower triangular 
IGTs in the reduction process. In Section |4.4[ this new understanding leads 
us to develop PREDUCTION, a new reduction algorithm which is more effi- 
cient and numerically stable than LAMBDA reduction and MREDUCTION. 



In Section |4T5| we present simulation results. Finally, in Section W7f\ we discuss 
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the implications of these results to the standard form of the OILS problem and 
to the LLL reduction algorithm. 

4.1 Impact of decorrelation on the search process 

As seen in Section 13.14 according to the literature, one of the two goals of 
the reduction process is to decorrelate the ambiguities as much as possible. 
Decorrelating the ambiguities as much as possible implies making the covari- 
ance matrix as diagonal as possible, i.e., making the absolute values of the 
off-diagonal entries of L as small as possible. In the following, we show that 
to solely make the absolute values of the off-diagonal entries of L as small as 
possible will have no impact on the search process. 

THEOREM 4.1.1. Given the OILS problem (J3TTJ) and the reduced OILS 
problem ( 13. 2\\ . If the transformation matrix Z is a product of lower triangular 
IGTs, then the search trees for problems ( 13. ip and (13. 2p are identical. 

Proof. Let the L T DL factorization of W& and Wz be 

W & = L T DL, W z = L T DL. 

As shown in Section I37T1 the OILS problems ( 13. ip and ( 13. 2ft can be written in 
the form (c.f. [3719]) 

n n 

/(aO = Z>i-*i)Vdi. /(*) = - hffii, (4- 1 ) 

3=1 3=1 

where 

n n 

Xj Xj ~\~ ^ lijij^i j %j Zj ~t~ ^ liji^Zi • (^'^) 

i=j+l i=j+l 

We first consider the case where Z is a single lower triangular IGT Z^j (k > j), 
which is applied to L from the right to make as small as possible (see 
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Section I3.1.ip . We have 

L = LZ kj = L - /ile^e] , 
where the modified entries in L are 

hj = hj - filth, t — k,...,n. (4.3) 

Let z = Z T kj x. Thus, 



Xi, if i ^ j, 

xj - x k fi, if i = j. 

From ffl~2]) and (jO|) . we know that 



(4.4) 



Zj = Xi, for i > j. (4.5) 

Lower triangular IGTs do not affect the D factor, meaning 

di = d h Vz. (4.6) 

We want to compare the enumerated points in the search process of problems 
( 13. ip and ( 13. 2D . The search process starts at level n and moves down to level 
1. When it moves down to level i, it chooses Xi = [xi] and zi = \_Zi\. From 
(14. 5 p and (14. 6p . if the chosen integer x^ is not valid, i.e., it does not satisfy 
bound ( I3.2ip at level i, then the chosen integer Z{ is also not valid. In this 
case, the search trees for problems (13. ip and (I3.2p will both move up to level 
2 + 1. Therefore, before we reach level j in the search process, we have 

^ = Xi, for % > j. 
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At level j, 



Xj Zj Xj ~\- ^ liji^Xi X{j Zj ^ ^ lij^i ^i) 

i=j+l i=j+l 



= X k /J + ^ l ij( X i ~ X i)- li i( Zi ~ Zj ) ( Usill § (S3D) 

i=j+l i=j+l 
n 

= x k fi + 2J (kj - kj)(%i - Xi) (using (SIS])) 

i=i+i 

n 

= x k fi + y^ y (kj - kj){xi - %i) (using (jl3D) 

i=k 
n 

= x k [i + ^2 rfikixi - %i) (using (@~3])) 

i=k 

n, 

= X k fX + [il kk {x k - X k ) + lll ik (Xi-Xi) 

i=k+l 

n 

= x k \i + /i[x fc + /^(xi - Xj) - S fc ] (since Z fefe = 1) 

i=fe+l 

= x kf j, (using g^D). (4.7) 

Since x k and fi are integers, Zj is an integer distance from Xj. This means that 
if integer Xj is chosen when we move down to level j in the search, then the 
chosen integer Zj is (see (14. 7ft ) 

Zj = Xj — x k fi. (4.8) 
From (14. 7p and (14. 8p . we obtain 

z j ~~ — Xj — Xj. (4-9) 
Using (JO]) and (J49]) in (Q, we get 

= Xj, for i < j. 

In other words, while Zj and Xj have different values, they have the same impact 
on the lower levels of the search process. Hence, the enumerated points in the 
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search process satisfy 

Zi = Xi, for % < j. 

For each level % in the search process, the enumerated points x and z satisfy 

zi = Xi, Vz ^ j, 

Zj = Xj — Xfc/i, for i = j. 



This shows that the search trees for problems (13.11) and (13.21) are identical. 
Consider the case where Z is a product of lower triangular IGTs, i.e., Z = 
Zi, . . . , Z n , used to make the absolute values of the other off-diagonal entries 
of L as small as possible. As shown, applying Z% to (13. ip will transform the 
ILS problem, but not modify the search tree. Applying Z 2 to this transformed 
ILS problem will also not modify the search tree, and so on. Thus, if Z is a 
product of lower triangular IGTs, the search trees for problems (13. ip and (13. 2 p 
are identical. □ 

Since the search trees are identical, lower triangular IGTs by themselves 
have no impact on the search process. Hence, it is not true that the search 
process is more efficient when the off-diagonal entries of L are as small as 
possible. 

We provide a 2 by 2 example. Let the covariance matrix Wj be 



Ws 



^ 11026 1050 ^ 



1050 100 



Its L and D factors are 



' 1 0^ 



10.5 1 



D 



'l ^ 



100 
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Let the real least-squares estimate be x = (5.38, 18.34) T . We can make |/ 2 i| 
as small as possible with the following IGT 



V 



i o 

-10 1 



/ 



The covariance matrix and its L factor become 



z 1 w^z 



^ 26 50 ^ 



In 



V 



50 100 



LZ 



' 1 0^ 



/ 



0.5 1 



the correlation coefficient p and the elongation of the search space 
e are used to quantify the correlation between the ambiguities. The correlation 
coefficent p between random variables s\ and S2 is defined as (see [3_3, p. 322]) 

The elongation of the search space e is given by square of the condition number 
of the covariance matrix (see Section I4.6p . For the original ambiguities x, 
we have p = 0.999 and e = 1.113 x 10 3 . For the transformed ambiguities 
z, we have p = 0.981 and e = 12.520. These measurements indicate that 
the transformed ambiguities are more decorrelated. The points (xi,X2) T and 
{zx,Z2) T encountered during the search process are shown in Table FT1 and 
14-2 1 where — indicates that no valid integer is found. In both cases, the 
first point encountered is valid, while the others points are invalid. The OILS 
solution is x = (2, 18) T . As expected, we observe that the lower triangular 
IGT did not reduce the number of points encountered in the search process. 



4.1.1 Implications to some reduction strategies 



In 



26], a united ambiguity decorrelation approach is proposed: unlike the 



usual pairwise decorrelation, all the ambiguities are decorrelated at once. This 
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Table 4-1: Search process without IGT 





Xo 


2 


18 




18 




19 




17 




20 







Table 4-2: Search process with IGT 



Zi 


Z2 


-178 


18 




18 




19 




17 




20 







allows for faster, but not maximum, ambiguity decorrelation. Their approach 
can be divided in two stages: (i) reordering the ambiguities (ii) decorrelating 
them. The reduction process can be written as M T P T WxPM, where P 
is a permutation matrix and M is a product of lower triangular IGTs. In 
this contribution, we have shown that stage (ii) will not improve the search 
process. This means that the search process would have been identical if the 
reduction strategy consisted of (i) only. The united ambiguity decorrelation 
approach can be iterated until no more decorrelation is possible. In this case, 
only the last decorrelation step can be removed. 

In the lazy transformation strategy of the MREDUCTION algorithm (see 
Section |3~3]) . when no more permutations occur in the reduction process, lower 
triangular IGTs are applied to the off-diagonal entries of L. Our current 
understanding shows that these IGTs are unnecessary; removing them reduces 
the computational cost of the reduction process, without affecting the search 
process. 
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4.2 Partial reduction 



In the literature, it is often conjectured that when the ambiguities get more 
decorrelated, the computational cost of the search process decreases. We have 
already shown that to solely decorrelate the ambiguities by applying lower 
triangular IGTs to the L factor of the L T DL factorization of the covariance 
matrix will not help the search process. However, as in LAMBDA reduction, 
lower triangular IGTs combined with permutations can significantly reduce 
the cost of the search process. This indicates that the accepted explanation 
is, to say the least, incomplete. We now provide a new explanation on the role 
of lower triangular IGTs in the reduction process. 

We claim that the computational cost of the search depends mainly on 
the D factor. The off-diagonal entries of L are only important when they 
affect D. In the reduction process, when we permute pair (k,k + 1), D is 
modified according to (13. 7ft . We strive for (13. 4ft . In order to make dk+i as 
small as possible, from (13.71) . we observe that |ife+i,fc| should be made as small 
as possible. An example would be helpful to show this. Let the L and D 
factors of a 2 by 2 covariance matrix W& be 





f i 


0^ 






\ 


L = 














^ 0.8 


1 J 






100 ) 



We have d>2 = 100 and d\ = 1. Let the real least-squares estimate be x = 
(13.5, 1.2) T . If we permute the two ambiguities without first applying an IGT, 
i.e., 

(o A 
V 1 
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using ( 13. 7p . we have d 2 = 65 and d\ = 1.54. The search process will be 
more efficient after this transformation because d 2 < d 2 allows more prun- 
ing to occur (see (13 .4f) ) . The integer pairs (zi,z 2 ) T encountered during the 
search are (2, 14) T , (-, 14) T , (-, 13) T and (-, -) T . The OILS solution is x = 
Z~ T (2, 14) T = (14, 2) T . However, we can make d 2 even smaller by applying a 
lower triangular IGT before the permutation, which means that 

(o 1 ^ 

Z = 

I 1 -v 

In this case, we have d 2 = 5 and d\ = 20. Now, three and not four integer pairs 
are encountered during the search, namely (2, 12) T , (— , 12) T and (— , — ) T . The 
OILS solution is x = Z~ (2, 12) T = (14, 2) T . This example illustrates how a 
lower triangular IGT, followed by a permutation, can prune more nodes from 
the search tree. 

It is useful to make |/fc+i,fc| as small as possible because of its effect on 
D. However, making as small as possible, where j > k + 1, will have 
no effect on D since ( 13.70 only involves h+i,k- Hence, even if is very 
large, making it smaller will not improve the search process. This means that 
making all the off-diagonal entries of L as close to as possible is unnecessary. 
It is only necessary to make |Zfc+i,k| as close to as possible before permuting 
pair (k, k + 1) in order to strive for (13.40 . We call this strategy a "minimal" 
reduction (MINREDUCTION) strategy. The MINREDUCTION algorithm is 
exactly like Algorithm 13.2. 1\ except that the first "for loop" is replaced with 
the following statement: [L, x, Z] = GAUSS(-L, k + 1, k, x, Z). 

Large off-diagonal entries in L indicate that the ambiguities are not decor- 
related as much as possible, which contradicts the claim that it is one of the 
goals of the reduction process. In the following, we provide a 3 by 3 example 
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which illustrates the issue. Let the covariance matrix and the real least-squares 
estimate of the ambiguity vector x be 



^ 2.8376 -0.0265 -0.8061 ^ 



\ 



-0.0265 0.7587 2.0602 
-0.8061 2.0602 5.7845 



f 26.6917 ^ 



x 



J 



\ 



64.1662 
42.5485 



/ 



Let ip(W x ) denote the sum of the absolute values of the correlation coeffi- 
cients of W x , which is used to quantify the entire correlation between the 
ambiguities. It is defined as follows 

n 

r(W. r ) = J2 Wj;.. l ) x /Wji.i)Wj.,.j) . 

i,j>i 

For the original ambiguities x, we have ip{W x ) = 1.2005. With LAMBDA 
reduction or MREDUCTION, the transformation matrix is 



/ 



V 



4 -2 1 
-43 19 -11 
16 -7 4 



\ 



/ 



The covariance matrix and the real least-squares estimate become 



/ 0.2282 0.0452 -0.0009 N 



Z T W,Z 



\ 



0.0452 0.1232 -0.0006 
-0.0009 -0.0006 0.0327 



Z'x 



\ 



-1971.6 
867.9 
-508.9 



\ 



We have ip{W z) = 0.2889, which indicates that the transformed ambiguities 
z are less correlated than the original ambiguities. With MINREDUCTION, 
the transformation matrix is 



/ 



Z = 



1 

1 -3 -11 
1 4 



\ 
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The covariance matrix and the real least-squares estimate become 



^ 0.7587 -0.2160 -0.1317 ^ 



W, 



z 1 w*z 



\ 



-0.2160 0.2518 0.0649 
-0.1317 0.0649 0.0327 



) 



( 



\ 



64.1662 
-149.9499 
-508.9418 



Now, we have ip{W z ) = 2.0449, which means that the transformed ambiguities 
are more correlated than the original ambiguities. 

Table 4-3: Search process with NOREDUCTION 





Z2 


Z3 


23 


64 


43 




64 


42 


27 


64 


42 




64 


42 






44 






41 









Table 4-4: Search process with LAMBDA reduction or MREDUCTION 



Zl 


Z 2 


Z3 


-1972 


868 


-509 




868 


-509 









Table 4-5: Search process with MINREDUCTION 



Zl 


z 2 


z 3 


64 


-150 


-509 




-150 


-509 









We refer to the reduction process which consists of only finding the L T DL 
factorization of W-j. as NOREDUCTION. The integer triples encountered dur- 
ing the search when NOREDUCTION, LAMBDA reduction and MINREDUC- 
TION are used are shown in Tables 14-3 [ 14-41 and 14-5 [ where — indicates that 
no valid integer is found. The ILS solution is x = (27, 64, 42) T . Observe 

42 



that MINREDUCTION causes the search to encounter the same number of 
integer points as LAMBDA reduction. Furthermore, NOREDUCTION causes 
the search to encounter four extra integer points than MINREDUCTION, al- 
though the ambiguities transformed by the latter are more correlated than 
the original ambiguities. This indicates that it is not true that the reduction 
process should decorrelate the ambiguities as much as possible in order for the 
search process to be more efficient. 

The significance of this result is threefold: 

• It indicates that contrary to common belief, the computational cost of 
the search is largely independent of the off-diagonal entries of L and of 
the correlation between the ambiguities. 

• It provides a different explanation on the role of lower triangular IGTs 
in the reduction process. 

• It leads to a more efficient reduction algorithm, see Section 14.41 



4.3 Geometric interpretation 



In Section I4.1[ we have shown that solely decorrelating the ambiguities will 
not improve the search process. We now illustrate this result geometrically. 
Let the real LS estimate of x and the covariance matrix be 



Xi 



X 



( 



V 



<J x 2 x 1 Of 2 



/ 



\ X2 / 

We assume that |cr£ i: £ 2 | > \& 2 Xl \ otherwise, no further decorrelation is possible. 
From f )3.20p . we observe that the ambiguity search space is centered at x. To 
decorrelate the two ambiguities, we use the following IGT 

\ 



/ 

V 



o 



'x 1 x 2 ^ X 



J 
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Note that the Z-transformation reduces of but does not affect cxj 2 in 



Therefore, as explained in 
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p. 365], the result of this transformation is to 
push the vertical tangents of the ambiguity search space, where the vertical 
axis is X2 and the horizontal axis x\ (see Fig. 14-1 p . Since |det(Z)| = 1, a 
Z-transformation does not change the area of the search space. 




Figure 4-1: Original and transformed search space. The transformed search 
space is less elongated. 



In Fig. 14-1 [ we observe that if the integer at level 2 is determined (vertical 
axis), the number of valid integers at level 1 (horizontal axis) is the same in 
the orginal and transformed search space. We prove this result as follows. In 
the search process, the lower bound 4 and the upper bound u k of the valid 
integers at a given level k are given by (13.221) and (13.231) . In the proof of 
Theorem I4.1.1[ we have shown that after applying an IGT Zjk for j > k, we 
have 

^ - %i)2 = fe"*') 2 , V z > k. 
di di 

We have also shown that Zk is an integer distance S of Xf.. These results imply 
that interval lift] for Xk and interval [Ik, Uk] for Zk satisfy (see ( 13.221) and 
IS!) 

[k,u k ] = [l k + 5,u k + 5]. 
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Hence in the search process, the number of valid integers Xk and the number 
of valid integers zj~ are the same. In other words, lower triangular IGTs by 



themselves do not reduce search halting 



In the literature (see, e.g. 



(see Sect.EZD. 



37]). it is commonly stated that the 



and 

elongated shape of the ambiguity search space (see (13.201) ). which is an hyper- 
ellipsoid, causes the search to be highly inefficient. It is then said that the role 
of the Z'-transformation is to make the search space less elongated. We now 
clarify this explanation. The elongation of the search space is defined to be the 
ratio of the major and minor principal axes. While a lower triangular IGT will 
not improve the search, it will, however, reduce the elongation of the search 
space (see Fig. 14-11 and the example given in Sect. 14. ip . This is explained by 
the fact that changes in the principal axes can occur without changes in the 
conditional variances. This shows that using the elongation of the search space 
to evaluate the effectiveness of the reduction can be misleading. See Section 
14.61 for an example and the relationship between the elongation of the search 
space and the condition number of the covariance matrix. 

4.4 A new reduction method 

This new understanding on the role of IGTs in the reduction process led us 
to a new algorithm: Partial Reduction (PREDUCTION). The motivation of 
PREDUCTION is to eliminate the unnecessary IGTs applied in LAMBDA 
reduction. Our results indicate that lower triangular IGTs have two roles in 
the reduction process. 

Efficiency for search. In Section 14.21 we showed that if pair (k, h + 1) will 

be permuted, then we first need to make 14+i.fcl as small as possible in 
order to improve the search efficiency. 
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Stability for reduction. When only making < 1/2, it is possible that 

large entries will appear in L(k + 2 : n, k). This effect may accumulate 
during the reduction process and introduce huge numbers, which can 
cause serious rounding errors. This problem does not occur if every time 
after we reduce |ifc+i,fe|, we also reduce |ife+2,fe|, • • • , \l n k\- 
The MINREDUCTION strategy mentioned in Section 14.21 did not include 
the IGTs that are necessary to ensure numerical stability. This means that 
on some ILS problems, MINREDUCTION yields a different and worse ILS 
solution than LAMBDA reduction. The PREDUCTION algorithm can be 
summarized as follows. It starts with column n — 1 of L. At column k, it 
computes the new value of lk+i,k if an IGT were to be applied. Then, using this 
new value, if d^+i > d k+1 holds (see (13.71) ). then permuting pair (k, k + 1) will 
not help strive for (13. 4ft , therefore it moves to column k — 1 without applying 
any IGT; otherwise it first applies IGTs to make \Uk\ < 1/2 for i = k + 1 :n, 
then it permutes pair (k, k + 1) and moves to column k + 1. In the LAMBDA 
reduction algorithm (see Sect. 13.21) . when a permutation occurs at column k, 
the algorithm restarts, i.e., it goes back to the initial position k = n — 1. From 
(13. 6p . no new permutation occurs in the last columns n — k — 1 of L. Hence, 
the new algorithm does not restart, but simply moves to column k + 1. In 
order to further reduce the computational costs, we use the symmetric pivoting 
strategy presented in Section 13.3.11 

We now present the complete PREDUCTION algorithm: 
Algorithm 4.4.1. (PREDUCTION). Given the covariance matrix and 
real-valued LS estimate x of x. This algorithm computes an integer unimod- 
ular matrix Z and the L T DL factorization W% = Z T WxZ = L T DL, where 
L and D are updated from the factors of the L T DL factorization of W±. This 
algorithm also computes z = Z T x, which overwrites x. 
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function: [Z, L, D, x] = PREDUCTION(VF 5; , x) 
Compute the L T DL factorization of W x 
with symmetric pivoting P T W X P = L T DL 
x = P T x 
Z = P 
k = n — 1 
kl = k 
while k > 

I = L(k + 1, k) - [L(k + 1, k)]L(k + 1, jfc + 1) 
D(k + 1, k + 1) = D(Jfe, fc) + / 2 D(A: + 1, jfc + 1) 
if D(k + 1, k + 1) < + 1, jfc + 1) 
if jfc < ifel 

for i = k + 1 : n 
II See Alg. I3XT1 
[L, x, Z] = GAUSS(L, i, k, x, Z) 

end 
end 

// See Alg. 13X21 

[L, D, x, Z] = PERMUTE(i, D, jfc, D(k + 1, k + 1), x, Z) 

kl — k 

if k < n — 1 
jfc = £; + 1 

end 
else 

k = k-1 
end 

end 
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Note that our final L might not be LLL-reduced since we do not ensure 
property (13. lip . Consider PREDUCTION without the symmetric pivoting 
strategy. Then, Theorem 14 . 1 . 1 1 implies that PREDUCTION will have the same 
impact on the search process as LAMBDA reduction. With the symmetric 
pivoting strategy, the initial ordering of the columns L is different than in 
LAMBDA reduction. For this reason, it is no longer true that the search 
process will be identical. Nevertheless, we do not expect significant differences 
in the computational cost, which is confirmed by simulations (see Sect. 14.51) . 
Unlike MREDUCTION, PREDUCTION ensures that we do not create large 
off-diagonal entries in L during the reduction process. This is necessary to 
avoid serious rounding errors. 

In the following, we give a 3 by 3 example to illustrate the differences 
between LAMBDA reduction, MREDUCTION and PREDUCTION. Let the 
covariance matrix and the real least-squares estimate of the ambiguity vector 
x be 



^ 1.3616 1.7318 0.9696 ^ 



V 



1.7318 2.5813 1.4713 
0.9696 1.4713 0.8694 



^ 27.6490 ^ 



x 



\ 



10.3038 
5.2883 



For this example, we get the same transformation matrix from LAMBDA 
reduction and MREDUCTION, which is 



-1 


1 





2 





1 


-2 


-1 


-2 
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The covariance matrix and the real least-squares estimate become 



^ 0.3454 -0.0714 0.0200 ^ 



W, 



z 1 w*z 



\ 



-0.0714 0.2918 0.0601 
0.0200 0.0601 0.1738 



z = Z T x 



' -17.6179 ^ 



22.3607 
-0.2727 



The L and D factors of the L T DL factorization of W \ are 

/ 1 n n\ / 



V 



1 

-0.2889 1 
0.1149 0.3459 1 



D 



\ 



0.3205 
0.2710 
0.1738 



\ 



/ 



The integer triples encountered during the search are shown in Table 14-6 \ 
where — indicates that no valid integer is found. The last full integer point 
found is z = (—18, 23, 0) T . The integer least-squares solution x for the original 
ambiguities is x = Z~ z = (28, 10, 5) T . Now, we compare the results with 

Table 4-6: Search process with LAMBDA reduction or MREDUCTION 



Zl 


Z2 


z 3 


-17 


22 







22 





-18 


23 












our PREDUCTION method. The transformation matrix is 

^01 ^ 

1 

1 -1 -2 



V 



The covariance matrix and the real least-squares estimate become 



^ 0.8694 0.1002 -0.2676 X 



Wi = z 1 w^z 



0.1002 0.2918 0.0601 
-0.2676 0.0601 0.1738 



Z T x 



( 5.2883 ^ 



22.3607 
-0.2727 
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The L and D factors of the L T DL factorization of W z are 

/ 1 n n\ / 



V 



1 

0.7111 1 
-1.5392 0.3459 1 



D 



J 



\ 



0.3205 
0.2710 
0.1738 



The integer triples encountered during the search are shown in Table 14-71 
The last full integer point found is z = (5, 23, 0) T . The integer least-squares 
solution x for the original ambiguities is x = Z~ T z = (28, 10, 5) T . Notice that 
in this example, the search process is identical whether the ^-transformation 
comes from LAMBDA reduction, MREDUCTION or PREDUCTION. 

Table 4-7: Search process with PREDUCTION 





Z2 




5 


22 







22 





5 


23 












4.5 Numerical simulations 

We implemented the PREDUCTION method given in Section 14.41 We did 
numerical simulations to compare its running time with LAMBDA reduction 
and MREDUCTION. All our computations were performed in MATLAB 7.9 
on a Pentium-4, 2.66 GHz machine with 501 MB memory running Ubuntu 
8.10. 

4.5.1 Setup 

We performed simulations for different cases. Cases 1-8 are test examples 



10j. With the exception of case 9, the real vector x was constructed 



given m 
as follows: 

x = 100 * randn(n, 1), (4.10) 
50 



where randn(n, 1) is a MATLAB built-in function to generate a vector of n 
random entries which are normally distributed. 

The first four cases are based on W& = L T DL where L is a unit lower 
triangular matrix with each 1^ (for i > j) being a random number generated 
by randn, and D is generated in four different ways: 

• Case 1: D — diag(cii), di = rand, where rand is a MATLAB built-in 
function to generate uniformly distributed random numbers in (0, 1). 

• Case 2: D = diagfn" 1 , (n - . . . , l" 1 ). 

• Case 3: D = diag(l _1 , 2" 1 , . . . ,n" 1 ). 

• Case 4: D = diag(200, 200, 200, 0.1, 0.1, . . . , 0.1). 
The last five cases are as follows: 

• Case 5: W A = UDU T , U is a random orthogonal matrix obtained 
by the QR factorization of a random matrix generated by randn(n, n), 
D = diag(dj), di = rand. 

• Case 6: Wx = UDU T , U is generated in the same way as in Case 5, 
di — , d n — 2% , other diagonal elements of D is randomly distributed 
between d 1 and d n , n is the dimension of W x . Thus the condition 
number of W& is 2t 

• Case 7: = A T A, A = randn(n,n). 

• Case 8: = UDU T , the dimension of is fixed to 20, U is 

k k 

generated in the same way as in Case 5, d\ = 2~?, d n = 2a, other 
diagonal elements of D are randomly distributed between d\ and d n , 
k = 5, 6, . . . , 20. Thus the range of the condition number of W& is from 
2 5 to 2 20 . 

• Case 9: We assume the linear model 

y = Ax + v, 
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where A = randn(2n, n), x = [100 * randn(n, 1)] and v ~ A/"(0, 0.05/). 
Then, we solve the following OILS problem 



min \\y — Ax\\ 2 



which we can rewrite in terms of (I3.ip ; see (11. 7p . 
Case 4 is motivated by the fact that the covariance matrix W x i n GNSS 
usually has a large gap between the third conditioned standard deviation and 



the forth one (see 
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Sect. 8.3.3]). The motivation for case 9 is that in 
typical GNSS applications, the variance of the noise vector v is small. For the 
reduction process, we took dimensions n = 5:40 and performed 40 runs for the 
all cases. For the search process, we took dimensions n = 5 : 30 and performed 
40 runs for all cases. The results about the average running time (in seconds) 
are given in Figs. 14-21 to 14-101 For each case, we give two plots, corresponding 
to the average reduction time and the average search time, respectively. Note 
that Z T W X Z = L T DL. Thus, W x = Z T L T DLZ 1 is a factorization of 
W x . We use the relative backward error to check the numerical stability of 
the factorization, which is 

\\W x -Z: T L^D c L c Z-% 
\\W X \\ 2 

where Z c , L c and D c are the computed values of Z, L and D. The results for 
the three reduction algorithms are displayed in Figs. 14-111 to 14-191 

4.5.2 Comparison of the reduction strategies 

From the simulation results, we observe that PREDUCTION improves the 
computational efficiency of the reduction stage for all cases. Usually, the 
improvement becomes more significant when the dimension n increases. For 
example, in Case 3, PREDUCTION has about the same running time as 
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MREDUCTION and LAMBDA reduction when n = 5, but is almost 10 times 
faster when n = 40. 

Below, we show that the MREDUCTION is not numerically stable. For 
this reason, we did not compare the effectiveness of MREDUCTION with the 
other reduction algorithms. With LAMBDA reduction and PREDUCTION, 
we obtain the same computed solution for the same OILS problem. 

In Section I4.4[ we showed that PREDUCTION without the symmetric 
pivoting strategy and LAMBDA reduction have exactly the same impact on 
the search process. With the symmetric pivoting strategy, PREDUCTION and 
LAMBDA reduction can give different D factors. The D factor which satisfies 
order (13. 4p better depends on the specific OILS problem. Nevertheless, Figs. 
14-21 to 14-91 show that there is no significant difference in the search process 
whether we use PREDUCTION or the LAMBDA reduction. 

In our simulations, we found that MREDUCTION sometimes causes the 
search to find a different and worse OILS solution than if LAMBDA reduction 
or PREDUCTION was used. For instance, this occured twice out of 180 runs 
for case 7 at n = 35. For these problems, the relative backward error of 
MREDUCTION was in the order of 10 2 and 10 12 , while the relative backward 
error of LAMBDA reduction and PREDUCTION was in the order of 10~ 14 . 
Observe that the relative backward error of MREDUCTION is particularly 
large for cases 2 and 7. This means that the MREDUCTION algorithm is not 
backward stable. The stability problem is caused by the lazy transformation 
strategy (see Section |3~3]) . The deferred IGTs in the reduction process can 
cause large off-diagonal entries in L to appear, which can lead to big rounding 
errors. Such an issue is avoided in PREDUCTION. For all the cases, the 
relative backward error in PREDUCTION is less than in LAMBDA reduction 
and in MREDUCTION. For some cases, the difference is of several orders 
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of magnitude. For example, for case 3 at n = 40, the relative backward 
error in LAMBDA reduction and MREDUCTION is around 1(T 12 , while it 
is 10~ 16 with PREDUCTION. This indicates that PREDUCTION is more 
computationally efficient and stable than both MREDUCTION and LAMBDA 
reduction. In some applications, the computational cost of the search process 
can be prohibitive. In order to reduce the search time, one might opt to find 
an approximate OILS solution. In these cases, the savings in the reduction 
time provided by PREDUCTION become particularly important. 
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Figure 4-2: Running time for Case 1 
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Figure 4-3: Running time for Case 2 
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Figure 4-4: Running time for Case 3 
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Figure 4-5: Running time for Case 4 
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Figure 4-6: Running time for Case 5 
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Figure 4-7: Running time for Case 6 
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Figure 4-8: Running time for Case 7 
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Figure 4-9: Running time for Case 8 
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Figure 4-10: Running time for Case 9 
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Figure 4—11: Relative backward error for Case 1 
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Figure 4-12: Relative backward error for Case 2 
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Figure 4-13: Relative backward error for Case 3 
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Figure 4-14: Relative backward error for Case 4 
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Figure 4-15: Relative backward error for Case 5 
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Figure 4-16: Relative backward error for Case 6 
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Figure 4-17: Relative backward error for Case 7 
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Figure 4-18: Relative backward error for Case 8 
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Figure 4-19: Relative backward error for Case 9 



4.6 Condition number criterion 



In some GNSS literature (see, e.g., [26|, [27J and |43j), it is believed that the 
goal of reduction process is to reduce the condition number of the covariance 



matrix W x . The 2-norm condition number of is defined as (see 20j, Sect. 
2.7]) 

«(W*) = ||W & ||||WT 1 || = cniW^/aniWz), 

where (Ji(W&) and cr n (W&) are the smallest and largest singular values of 
Wx- Geometrically, the condition number corresponds to the square of the 



ratio of the major and minor axes of the search ellipsoid (see 33|, Sect. 6.4]). 
In other words, the condition number measures the elongation of the search el- 
lipsoid. As seen in Section I4.3[ decorrelating the ambiguities makes the search 
ellipsoid less elongated. Thus, a lower condition number of the covariance 
matrix indicates that the ambiguities are more decorrelated. However, our 
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contribution has shown that the goal of the reduction process is not to decor- 
relate the ambiguities as much as possible. We now provide an example which 
shows that the condition number criterion can be misleading. Let 



W 6 



L T DL 



1 1000.5 
1 



\ 



/ 



4 
0.05 



\ 



/ 



V 



1 
1000.5 1 



/ 



The condition number of is 1.2527 x 10 10 . If we apply a lower tri- 
angular IGT to decorrelate the 1st and 2nd ambiguity, our new covariance 
matrix is 



L T DL 



( 

V 



\ 



4 
0.05 



\ 



\ 



1 
0.5 1 



1 0.5 

1 J \ U U.UO J \ U.O LI 

The condition number of W z is 80.5071. This shows that to solely decor- 
relate the ambiguities can drastically lower the condition number of the covari- 
ance matrix, yet it will not yield any improvement towards the search process 
as shown in Section I4TT1 If we permute the ambiguities, the L T DL factorization 
becomes 



P T W Z P = L T DL 



1 0.0062 
1 



\ 



/ 



/ 

V 



0.0498 
4.0125 



\ 



/ 



1 
0.0062 1 



/ 



The condition number of P T W %P is still 80.5071. Yet, numerical simu- 
lations indicate that the search process is faster with for randomly chosen 
x. This is explained by the fact that order f)3.4p is satisfied by and not 
P T W zP. For this reason, order H3.4[) is a better criterion to evaluate the 
reduction process than the condition number. 
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4.7 Implications to the standard OILS form 



We now translate our result on the role of lower triangular IGTs in the 
quadratic OILS form (13. ip to the standard OILS form (12. ip . Note that in 
the standard form, it is upper triangular IGTs that make the absolute values 
of the off-diagonal entries of R as small as possible. Section 14.11 showed that 
making R as close to diagonal as possible will not help the search process. We 
know that the purpose of permutations is to strive for ru -C . . . <C r nn . This 
contribution shows that it is also the purpose of upper triangular IGTs. In 



the reduction process, if rk-x,k-x > \lr\_ x k + r\ kl then permuting columns k 



and k — 1 will decrease Tk-i^-i and increase r^k (see Section [2.2.2j) . The pur- 
pose of an upper trianguh IGT is two-fold. First, applying IGT Zk-\ y k makes 
l r /c-i,fc| as small as possible, which increases the likelihood of a column per- 
mutation. Second, if a permutation occurs, from (12.181) and (12.191) . a smaller 
l r fe-i,fc| ensures a greater decrease in rk-i,k-i an d a greater increase in r&&. 

4.7.1 A partial LLL reduction 

Our result indicates (wrongly) that IGTs have to be applied only on the su- 
perdiagonal entries of R. Such an approach can be numerically unstable since 
it can produce large off-diagonal entries in R relative to the diagonal entry 
(see [23|). Some other IGTs are needed to bound the off-diagonal entries of 
R. The chosen approach is as follows. If columns k and k — 1 are going to be 
permuted, we first apply IGTs to make the absolute values of rk-i t k, • • • , rut as 
small as possible; otherwise, we do not apply any IGT. This way, the number 
of IGTs is minimized. Based on this idea, we present a partial LLL (PLLL) 
reduction algorithm. 

Algorithm 4.7.1. (PLLL Reduction). Given generator matrix A e K mxn 
and the input vector y G M m . The algorithm returns the reduced upper 
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triangular matrix R £ M rixri , the unimodular matrix Z e Z nxn , and the 
vector y £ IR n . 

function: [Ji, Z, y] = PLLL(A, y) 



Compute the sorted QR decomposition of A (see 
and set y = Wjy 
Z = I 
k = 2 

while k < n 

r t = r k -i,k - L^fc-i.Jfe/^fe— i,fc-il x r fc _i, fc _i 



if rfc-i^-i > A/r t 2 + 



kk 



for i — k — 1: — 1:1 

Apply IGT Z ife to .i.e., = i?Z, ;fc 
Update Z, i.e., Z = ZZjfc 

end 

Interchange columns k and k — 1 of R and Z 

Transform i? to an upper triangular matrix by a Givens rotation 

Apply the same Givens rotation to y 

if k > 2 
k = k - 1 

end 
else 

fc = A; + 1 
end 



end 



We point out that our final R might not be LLL-reduced since we do not 
ensure property (12. 15ft . 
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4.7.2 Ling's and Howgrave-Graham's effective LLL reduction 

The initial version of this thesis was submitted before finding the effective LLL 
reduction presented in 13j. We now point out the similarities and differences 
with the results derived in this chapter. Firstly, they show that IGTs do not 
affect the Babai point, while we show that IGTs do not affect the entire search 
process. Secondly, their modified LLL reduction algorithm uses Gram-Schmidt 
orthogonalization. Our modified LLL reduction algorithm uses Householder 
reflections and Givens rotation. The latter is more stable. Thirdly, our al- 
gorithm does not apply an IGT if no column permutation is needed as it is 
unnecessary. Finally, they do not take numerical stability into account, while 
we do. Specifically, our reduction algorithm uses extra IGTs to prevent the 
serious rounding errors that can occur due to the increase of the off-diagonal 
elements. 
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CHAPTER 5 

Ellipsoid-Constrained Integer Least Squares Problems 



In some applications, one wants to solve 

mm\\y - Ax\\ 2 2 , £ = {x e Z n : \\Ax\\ 2 2 < a 2 }, (5.1) 

where y G W 1 and A G M mxn has full column rank. We refer to J5.ip as 
an ellipsoid-constrained i nteg er least squares (EILS) problem. In [15], the 



V-BLAST reduction (see 
EILS problem. In 



19]) was proposed for the reduction process of the 
9j, it was shown that the LLL reduction makes the search 
process more efficient than V-BLAST. It was also noted that for large noise 
in the linear model (see (ll.ip ). the search process becomes extremely time- 
consuming both with V-BLAST and the LLL reduction. In Section 15.11 we 
show how to modify the search process given the ellipsoidal constraint based 



on the work in 



9j. In Section |572| we give a new reduction strategy to handle 



the large noise case, which unlike the LLL reduction and V-BLAST, uses all 



the available information. Finally in Section 15. 3[ we present simulation results 
that indicate that our new reduction algorithm is much more effective than 
the existing algorithms for large noise. 

5.1 Search process 

Suppose that after the reduction stage, the original EILS problem (15.11) is 
transformed to the following reduced EILS problem 



min \\y - Rz\\l, £ = {z e Z n : \\Rz\\ 2 2 < a 2 }. (5.2) 

Z: f 
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Without loss of generality, we assume that the diagonal entries of R are posi- 
tive. Assume that the solution of (12.21) satisfies the bound 

\\y-Rz\\ 2 2 <(3 2 . 

Then, as in the OILS problem (12.21) . we have the following inequalities (see 
Sect. EI]) 

level n : r 2 nn (z n - c n f < [3 2 , (5.3) 



n 

level k : rl k (z k - c k f < f3 2 - £ r 2 u {z t - c,) 2 (5.4) 

i=k+l 



level 1 : r 2 u ( Zl - Cl ) 2 < f3 2 - $> 2 (^ - q) 2 , (5.5) 

where c k , for A; = n : 1, are defined in (12.41) . In Section 12.11 we presented 
a search process for the OILS problem based on these inequalities. For the 
EILS problem (I5.2p . the search also needs to take the constraint ellipsoid into 
account. In 9(, the search given in Section [2TT1 was modified in order to ensure 
that the enumerated integer points satisfy the constraint ellipsoid. Here, we 
show how to compute the bounds of the constraint. 

The constraint ellipsoid ||i?2|| 2 < a 2 can be written as 

n n 

^2(r kk z k + ^ r kjZj) 2 < a 2 (5.6) 

k+l j=k+l 

Define 

n 

b n = 0, b k = ^ r kjZj, k = n - 1 : -1 : 1. (5.7) 

j=k+i 

n n 

s n = a 2 , s fe _i = « 2 -^(rii^+ ^2 r ikZjf = s k -(r kk z k +b k ), k = n : -1 : 2. 



i=k j=i+l 



(5.8) 
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With flEZD and (EHJ, we rewrite (ESD as 



(^fcfc^fc + b k ) 2 < s k , k = n:-l:l. 



Therefore, at level k in the search process, z k is constrained to the interval 



rkk 



Uk 



y/Sk-h 



rkk 



k = n : — 1 : 1. 

(5.9) 

At each level in the search process, we compute Ik and Uk with (15. 9p . If 4 > 
then no valid integer exists at level A; and we move up to level k + 1. In the 
following, the Schnorr-Euchner search algorithm is modified to ensure that Zk 
is constrained to [lk,Uk\; see {9]. 

Algorithm 5.1.1. (SEARCH-EILS) Given nonsingular upper triangular ma- 
trix R G M. nxn with positive diagonal entries, the vector y e M n , the initial 
search ellipsoid bound (3 and the constraint ellipsoid bound a. The search 
algorithm finds the solution z € Z n to the EILS problem (15. 2p . 
function: z = SEARCH_EILS(i?, y, (3, a) 

1. (Initialization) Set k = n,bk = and = a 2 

2. Set Iboundk = and uboundk = 0. 



Compute Zfc = 
if M fc < l k 

Go to Step 4 
end 

if u k = h 

Set Iboundk - 
end 



ffcfc 



i'kk 



1 and uboundk 



Compute c fc = - b k )/r kk - Set z fc = |c fc "|, 
if z fc < 4 

z k = h, set Iboundk = 1 and A fc = 1 
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else if Zk > Uk 

Zk = Uk, set uboundk = 1 and A k = — 1 
else //no bound of the constraint is reached 

Set A k = sgn(c fe - z k ) 
end 

3. (Main step) 

go to Step 4 
else if k > 1 

Compute = Yl"=k r k-i,j z j, s k-i = s k - (r kk z k + b k ) 2 

Set k — k — 1, go to Step 2 
else / / case k = 1 

go to Step 5 
end 

4. (Invalid point) 
if k — n 

terminate 
else 

k — k + 1, go to Step 6 
end 

5. (Found valid point) 

Set z = z, p = ELi r L(4 - c k f 
k — k + 1, go to Step 6 

6. (Enumeration at level k) 

if uboundk = 1 and Iboundk = 1 

Go to Step 4 // no integer is available at this level 
end 
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Set z k = z k + A k 
if z k = l k 

Set Iboundk = 1, Compute Afc 
else if Zfc = Ufe 

Set uboundk = 1, Compute A fc 
else if Wounds = 1 

A fc = l 
else if uboundk = 1 

Afe = -1 
else 

Compute Afe = -A fe - sgn(A fc ) 
end 

Go to Step 3. 
5.2 Reduction 

As in the OILS problem, we can apply IGTs and permutations (see Section 12^21) 
in the reduction process of the EILS problem. With fj2TTU]) . (I2TTT) and (EHSjl . 
the original EILS problem (15. ip is transformed to the reduced EILS problem 
( 15. 2ft . Our new reduction for the EILS problem is based on a reduction strategy 
first applied to the box-constrained integer least squares problem. 

5.2.1 A constraint reduction strategy 

In several applications, one wants to solve 

mm\\y - Ax\\ 2 2 , B = {x e Z n : I < x < u, I G Z n , u G Z n }. (5.10) 
We refer to ( I5.10p as a box-constrained integer least squares (BILS) problem. 



= -Afc - sgn(Afc) 
= -Afc - sgn(Afc) 
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The transformations applied on A during the reduction can be described 
as a QR factorization of A with column pivoting: 



Q T AP 



R 




or AP = QTR, (5.11) 



where Q = [Qi,Q 2 ] £ M. mxm is orthogonal, R G lR n is nonsingular upper 
triangular and P G Z™ xri is a permutation matrix. This is a special case of 



the QRZ factorization presented in Section 12. 2\ where we have a permutation 
matrix P instead of a general unimodular matrix Z. The reason is that a 
general unimodular matrix will make the box constraint B very difficult to 
handle in the search. Hence, the LLL reduction is usually not used to solve 
the BILS problem ( 15.101) in the literature. With (15. lip , we have 



\\y 

Let 



- Ax\\l = WQly - RP T x\\l + \\Qly\\l (5.12) 



y = Q^y, z = P T x, I = P T l, u = P T u. (5.13) 
Then, from ( 15. 12ft and ( 15. 13ft . we see that ( 15.1 Oft is equivalent to 

min||j/ — -R«||l, B = {z G Z n : 7 < z < u,l G U\u G II 1 }. (5.14) 

If z is the solution of the transformed BILS problem (15.141) . then x = Pz is 
the solution of the original BILS problem (I5.10p . 

Most reduction algorithms are solely based on A. In 
that using the information of y and the constraint in the reduction process can 
make the search process much more efficient. We call their reduction strategy 
a "constraint" reduction strategy. Here, we describe its main idea. In the 



, it was shown 
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search process at level i, we have 

n 

rl(zi - Q) 2 < - r lk( z k - c k )\ (5.15) 

k=i+l 

where c, is determined when z i+ i, . . . , z n are fixed (see (12.41) ). If we can reduce 
the search range of Zi for % = n, n — 1, . . . , 1, then the search will be more 
efficient. Notice that this can be achieved if 

• The right-hand side of (15.151) is as small as possible, which means that 
each r\ k {z k — c^) 2 is as large as possible. 

• ru is as large as possible. 

The constraint reduction strategy looks for a permutation of A such that 
\fkk{ z k — Ck)\ is as large as possible for k = n, n — 1, . . . , 1, where we also take 
into account that r k k should be as large as possible. The algorithm determines 
the columns of the permuted A from right to left. To determine the fcth 
column, it chooses from the remaining k columns the one that maximizes 
\ r kk{ z k ~ Ck)\- Now, one question that arises is how to choose z k . The natural 
approach is to set z k to be the nearest integer to c k in \lk,Uk}. This can yield 
the following problem. If Zk is very close to Ck, then \r k k(zk — Ck) \ is small even 
though Tkk is large. Since ru...r nn is constant (note that det 1//2 (A T A) = 
det(-R) = ru . . . r nn ), we might end up with a large ru for small index i and 
a small ru for large index i. This does not comply with our requirement; see 
the first sentence of the paragraph. On the other hand, if we choose Zk to 
be the second nearest integer to Ck in [lk,itk\, then \zk — Ck\ is always larger 
than 0.5. Thus, if is large, then \rkk(zk — Ck)\ is also large. Hence, the 
previous problem is avoided. Simulations in [8| indicate that choosing Zk to be 
the second nearest integer to in is more effective than other possible 

choices of Zk- 
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5.2.2 The ellipsoidal constraint 

We want a reduction for the EILS problem which uses all the available infor- 
mation in order to improve the search. The natural approach would be to use 
the constraint reduction strategy with the ellipsoidal constraint. Our exper- 
iments show that such an approach is inefficient. Here, we explain why the 
constraint reduction strategy is effective for the BILS problem, but ineffective 
for the EILS problem. In the reduction process, we permute A such that its 
kth column maximizes \rkk{zk — Ck)\, where Ck depends on the chosen values of 
. . . , z n . In the EILS problem and unlike the BILS problem, the constraint 
depends on z (see (I5.2p ). In the search process at level k, Zk can take values 
in the interval [Ik, where Ik and Uk depend on Zk+i, . . . ,z n . As Zk+i, ■ ■ ■ , z n 
take on different values in the search, it is possible that the fcth column of A 
no longer maximizes \rkk(zk — Cfc)|. Numerical experiments indicate that if we 
have a box-constraint, it is likely that \rkk{zk~ Ck) \ remains large, which means 
that the search process remains efficient. If we have an ellipsoidal constraint, 
it is much less likely that \rkk(zk — Ck) \ remains large since in addition to Ck, 
the constraints Ik and Uk also change as Zk+i, ■ ■ ■ , z n change. In other words, 
the extra uncertainty in the EILS problem makes it more difficult to determine 
which column maximizes \rkk{zk — Qc)|- 

To overcome this difficulty, we construct the smallest hyper-rectangle 
which includes the constraint ellipsoid. The edges of the hyper-rectangle are 
parallel to the z-coordinate system. We suggest to use this new box-constraint 
instead of the constraint ellipsoid in the reduction. In the constraint reduction 
strategy, IGTs are not used since they make the box-constraint too difficult 
to handle in the search. This difficulty does not occur with the ellipsoidal 
constraint as shown in Section 15.11 Hence, we can modify the constraint re- 
duction strategy by introducing IGTs in the reduction stage for the EILS 
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problem. Note that the shape of the constraint ellipsoid changes after IGTs 
are applied, which means that the box-constraint needs to be recomputed. 
While the box-constraint is less precise than the constraint ellipsoid, it has 
the advantage of being insensitive to the chosen values of z in the reduction. 
This means that it is now more likely that \r k k{zk — c&)| remains large in the 
search process. 

5.2.3 Computing the box-constraint 

In 9] , it was shown how the smallest box-constraint [Z , u] that includes the 
constraint ellipsoid can be efficiently computed. For k = 1 : n, we want to 
determine 

Ik = Tmin eTz], given II Rz II 2 < a, (5.16) 

261" 

Uk = [maxe^zj, given ||-Rz|| 2 < a. (5-17) 
We first solve for u k . Let p = Rz. Substituting in (15. 17ft . we obtain 

u k = Lmaxe^i2 _1 pJ, given ||p|| 2 < a. (5.18) 
p 



Using the Cauchy-Schwarz inequality (see [20, p. 53]), 



e T k R^p < ||ir T e fc || 2 ||p|| 2 < ||^- T e fe || 2 «. (5.19) 

The inequalities become equalities if and only if p and R~ 7 'e& are linearly 
dependent, and ||p|| 2 = a, i.e., p = aR~ T e k /\\R~ T ek\\2- Substituting p in 
(15.18j) . we get u k = La||-R _T efc|| 2 J . Note that m&x ze ^n e k z = — mm ze s.n e k z 
implies that If. = \— a||i?~ T efc|| 2 ] . To efficiently compute R 7 e k) we solve for 
q in the lower triangular system R T q = e k . The algorithm is summarized as 
follows. 

Algorithm 5.2.1. (BOX). Given the nonsingular upper triangular R e M nxn , 
the constraint ellipsoid bound a and an integer k, where k — 1 : n. The 



algorithm computes the interval of the hyper-rectangle \l,u] which 

includes the constraint ellipsoid. 

function: [lk,Uk\ = BOX(R,a,k) 

Solve R T q = for q by forward substitution 

Compute Uk = |_ a iM|2j an d h = \— ®\\q\\2~\ 

5.2.4 A new reduction algorithm 

Our new algorithm for the EILS problem merges the ideas of the constraint 
reduction strategy, which uses all the available information, and of the LLL 
reduction, which applies IGTs to strive for rn < ... < r nn . We call our 
algorithm Constrained LLL (CLLL) reduction. In the CLLL reduction, we use 
constraint information and IGTs to strive for |rn(zi — ci)| < . . . < \r nn (z n — 

C n )\- 

We now describe our reduction algorithm. The CLLL reduction starts 
by finding the QR decomposition of A by Householder transformations, then 
computes y and works with R from left to right. At the /cth column of R, the 
algorithm applies IGTs to ensure that {r^l < \ra for i — k — 1 : — 1 : 1. Then, 
it computes the constraints of Zk (see Section I5.2.3P and approximates 

Ck- The reason that Ck needs to be approximated is that at the kth column, 
Zk+i, ■ ■ ■ , z n are not yet determined (see f l2.4p ). The algorithm approximates 
(12 .4p by Ck = Vk/ r kk for k = n: —1 : 1. As in the constraint reduction strategy 
(see Sect. 15 .2 . 1 j) . it sets Zk to be the second nearest integer to Ck in \lk,Uk\- 
Then, if permuting columns k — 1 and k maximizes \rkk{zk — c~k)\, it does so, 
applies a Givens rotation to R from the left to bring R back to an upper 
triangular form, simultaneously applies the same Givens rotation to y, and 
moves down to column k — 1; otherwise it moves up to column k + 1. We 
present our implementation of the CLLL reduction. 
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Algorithm 5.2.2. (CLLL REDUCTION). Given the generator matrix A G 
R mxn , the input vector y G M m and the constraint ellipsoid bound a. The 
algorithm returns the reduced upper triangular matrix R G R nxra , the uni- 
modular matrix Z G Z nx ™, and the vector y G M n . 

function: [R, Z, y] = CLLL(A, y, a) 

Compute the QR decomposition of A and set y = Q\y 

Z :=I n 

k = 2 

while k < n 

for i — k — 1 : — 1 : 1 

Apply the IGT Z ik to R, i.e, R := RZ ik 

Update Z, i.e, Z := ZZ ik 
end 

R' := R 

Interchange columns k — 1 and k of R' and transform R to 
an upper triangular matrix by a Givens rotation, G 
y' = Gy 

1 1 Compute the box constraint of R for z k 
[l k ,u k ] =BOX(R,a,k) 

1 1 Compute the box constraint of R' for z' k 
[l' k ,u' k ]=BOX(R',a,k) 
II Compute \r kk (z k - c k )\ 
Cfc := y k /r kk 

el := y'k/ r 'kk 

Set z k to be the second nearest integer to c k on \l k ,u k } 
Set z' k to be the second nearest integer to c' k on [4,^] 

if Kfc(4-c5b)l > kfcfc(^-cfc)| 
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R:= R 

y ■= y' 

if k > 2 

k = k-l 
end 
else 

k = k + 1 

end 

end 

Note that the CLLL reduction algorithm moves from the left to the right 
of R, which explains why Zk+i, ■ ■ ■ , z n are not determined at the fcth column. 
It is straightforward to modify the LLL reduction algorithm for it to move 
from right to left, which allows us to determine Ck exactly. Surprisingly, sim- 
ulation results indicate that such an approach is less effective than the CLLL 
reduction. Further investigation is required to understand why \rkk{zk — c&)| 
is a better criterion than \rkk(zk — c k)\ to determine the permutation of A. 

5.3 Numerical simulations 

In this section, we implemented the CLLL algorithm given in Section 15.2.41 
We did numerical simulations to compare its effectiveness with the LLL re- 
duction. Algorithm 15.1.11 is used for the search process. All our simulations 
were performed in MATLAB 7.9 on a Pentium-4, 2.66 GHz machine with 501 
MB memory running Ubuntu 8.10. 
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5.3.1 Setup 

We took A to be n x n matrices drawn from i.i.d. zero-mean, unit variance 
Gaussian distribution. We construct y as follows 



y = Ax + v, 

where the noise vector v ~ A/"(0, cr 2 /). To generate x, we randomly pick an 
integer point inside some hyper-rectangle. Then, we set a = ||Aa;|| 2 . In Figs. 
15-11 to [5^5] we display the average CPU search time in seconds for a = 0.5 to 
a = 10. We took dimensions n = 5 : 30 for Figs. 15-H and I5-2"| n = 5 : 15 for Fig. 
!5-3[ n = 5 : 10 for Fig. 15-4 [ and performed 20 runs for each case. In Fig. !5-5[ 
we took n = 4 : 8 and performed 5 runs. The reason that the dimensions of 
the experiments get smaller as a gets larger is that the search process with the 
LLL reduction becomes extremely time-consuming. For instance, for typical 
problems with dimension n = 15 and a = 2, the search time with the LLL 
reduction is more than 10 3 s, while the search time with the CLLL reduction is 
around Is. In Fig. 15-61 we compare the Babai integer points corresponding to 
the LLL reduction and the CLLL reduction. We give the ratio between /3 at 
the Babai integer point and (3 at the EILS solution with different noise, where 
(3 = (3(z) = \\y — RzWl- We present the results for dimension n = 5 with 20 
runs. We obtain similar results for other dimensions. As shown in 9] , the cost 
of the search time dominates the cost of the whole algorithm when o > 0.5. 
For this reason, the figures do not take into account the reduction time, which 
is negligible. 

5.3.2 Comparison of the reduction strategies 

From the simulations, we see that the most effective reduction algorithm de- 
pends on the noise size. In Figs. 15-11 and !5-2[ we see that when the noise is 
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Figure 5—1: Average search time versus dimension, a = 0.5. 
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Figure 5-4: Average search time versus dimension, a = 4. 
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small, i.e, a < 1, the LLL reduction is more effective than the CLLL reduc- 
tion. When the noise gets larger, the CLLL reduction becomes much more 
effective than the LLL reduction. Notice that the improvement of the CLLL 
reduction over the LLL reduction becomes more significant with larger noise. 
For example, when n = 7, the CLLL reduction is slightly more effective than 
the LLL reduction when a = 2, but close to 1000 times more effective when 
a = 10. We observe that with the LLL reduction, the search time becomes 
more and more prohibitive as the noise gets larger. The CLLL reduction pro- 
vides considerable savings in the search time. For example when a = 2 and 
n = 8, the search time with the LLL reduction is more than 100s, while it 
is about Is with the CLLL reduction. When a = 4 and n = 10, the search 
time with the LLL reduction is more than 10 4 s, while it is about 10s with the 
CLLL reduction. 

We now explain why the LLL reduction is preferable over the CLLL re- 
duction in Figs. 15-11 and 15-21 When the noise is small, it is likely that the OILS 
solution is close the ellipsoidal constraint. As can be seen in Fig. 15-61 in these 
situations the Babai integer point found with the LLL reduction is usually 
very close to the EILS solution. With the CLLL reduction, the goal is to make 
\fkk{zk — c k )\ as large as possible for k — n, . . . , 1. While | rifc=i r kk\ is constant 
through the reduction process, | YYk=i r kk(z k — Cfc)| is n °t- Making \r kk (z k — c k )\ 
large for some k will not make \rjj{zj — Cj)\ smaller for j < k. It is possible 
that the CLLL reduction permutes A such that X]fc=i r L(- 2 fc — Ck) 2 ends up 
to be a large number. While such an ordering allows to prune more points in 
the search process, it also means that the Babai integer point found with the 
CLLL reduction is usually worse than the one found with the LLL reduction. 
A better Babai integer point makes the intersection between the search ellip- 
soid and the constraint ellipsoid smaller. This implies that the search process 
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with the CLLL reduction must find many points before it reaches the EILS 
solution, while with the LLL reduction, only very few points are found before 
the EILS solution. For large noise (see Figs. 15-31 to HP5]) . it is no longer true 
that the Babai integer point found with the LLL reduction is very close to the 
EILS solution. It is here where the extra search pruning that the CLLL re- 
duction provides becomes advantageous. Thus, the LLL reduction is the most 
effective reduction algorithm for small noise, which includes many communi- 
cations applications; while the CLLL reduction is by far the most effective 
reduction algorithm for large noise, which includes the applications where the 
linear model (see fll.ip ) is not assumed. 
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CHAPTER 6 
Summary and future work 



This thesis was concerned with solving the ordinary integer least squares 
(OILS) problem 

min \\y — Axil 2 ., 

where y G M n and A G ]j mxn h as f u n column rank. In the GNSS literature, 
one needs to solve the following quadratic form of the OILS problem 

minfcc — x) T W^. l (x — x), 

where x G W 1 is the real-valued least squares (LS) estimate of the double 
differenced integer ambiguity vector x G Z n , and W x G M. nxn is its covariance 
matrix, which is symmetric positive definite. 

There are two steps in solving an OILS problem: reduction and search. 
The main focus of this thesis was on the reduction step. 

In Chapter 01 we have shown that there are two misconceptions about the 
reduction in the literature. The first is that the reduction should decorrelate 
the ambiguities as much as possible. We have proved that this is incorrect: 
only some ambiguities should be decorrelated as much as possible. Our new 
understanding on the role of IGTs in the reduction process led to the PRE- 
DUCTION algorithm, a more computationally efficient and stable reduction 
algorithm than both LAMBDA reduction and MREDUCTION. The second 
misconception is that the reduction process should reduce the condition num- 
ber of the covariance matrix. We gave examples which demonstrate that the 
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condition number is an ineffective criterion to evaluate the reduction. Finally, 
we translated our result from the quadratic OILS form to the standard OILS 
form. Our new understanding on the role of IGTs in the LLL reduction algo- 
rithm led to the more efficient PLLL reduction algorithm. 

In Chapter O we discussed how to solve the ellipsoid-constrained integer 
least squares (EILS) problem 

min||y-i4x|||, £ = {x e Z n : \\Ax\\l < a 2 }, (6.1) 

where |/6l" and A G M mxn has full column rank. With the existing reduc- 
tion algorithms for the EILS problem, the search process is extremely time- 
consuming for large noise. We proposed a new reduction algorithm which, 
unlike existing algorithms, uses all the available information: the generator 
matrix, the input vector and the ellipsoidal constraint. Simulation results in- 
dicate that the new algorithm greatly decreases the computational cost of the 
search process for large noise. 

In the future, we would like to investigate the following problems: 
• Box-constrained integer least squares (BILS) problems often arise in 
communications applications. In the literature of BILS problems, IGTs 
are usually not applied in the reduction phase since they make the box 
constraint too difficult to handle in the search phase. We would like to 
see if the modified box constraint can be efficiently approximated by a 
larger constraint, which can be used easily in the search. The search pro- 
cess would have to be modified to ensure that the search ellipsoid is only 
updated if the integer point found satisfies the initial box-constraint. 
While the larger constraint makes the search less efficient, the reduction 
becomes more effective due to the IGTs. 
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• For the BILS and EILS problem, it was shown that using all the avail- 
able information in the reduction phase can be very effective. For the 
OILS problem, the LLL reduction is solely based on the generator matrix 
A. We would like to determine if using y can lead to a more effective 
reduction algorithm. 
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